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DESIGN OF SOME NUMERICAL GENERAL CIRCULATION EXPERIMENTS 


YALE MINTz!.2,3 
The Hebrew University of Jerusalem 


ABSTRACT 


The equations of motion are developed for a two-parameter model of the atmosphere with 
the friction and heating terms retained. An empirical surface stress coefficient is found 
by the diagnostic use of the equations. 


A heating function is specified in which the diabatic heating of the air depends on the 
earth—air temperature difference. It is shown that this heating function can be used in the 
equations to study the mean zonal and mean meridional circulations of the atmosphere 
on a land-covered planet and on an ocean-covered planet in which the initial ocean 
temperature is a function of latitude; and to study the complete general circulation 
of the atmosphere on an ocean-covered planet in which the initial ocean temperature 
is a function of latitude and longitude and on a planet covered partly by ocean and 
partly by land. 


INTRODUCTION 


The availability of high-speed large-capacity electronic digital computers makes it 
possible to study the general circulation of the atmosphere by numerical integration 
of the equations of motion. 


In the first section of this paper the equations of motion are developed for a two- 
parameter model of the atmosphere with the friction and heating terms retained. 
In this development the character of the vertical wind shear is not assumed but is 
derived from the given static stability of the air. In addition, the orographic effect 
is not prescribed in the usual way, but at each level the forced vertical motion is set 
equal to the product of the horizontal velocity and the slope of an orographically 
determined stream-surface in the direction of the velocity. 


1 On leave of absence from the University of California, Los Angeles. 
2 U.C.L.A. Department of Meteorology, Contributions to Meteorology, No. 39. 


3 The research reported in this paper was sponsored in part by the Geophysics Research Direc- 
torate of the Air Force Cambridge Research Center, Air Research and Development Command, 
United States Air Force, through the European Office, ARDC. 
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In the second section the mean surface stress coefficient is found, in an empirical 
way, by using the equations of motion diagnostically. Integrated over a polar cap 
this is the method of the zonal momentum budget, and it gives a ratio of surface 
stress to surface wind of 0.32 x 10-4t m~?s~! (0.32 dyne cm~? per m s~! geostrophic 
surface wind). 

In the third section it is shown that the complex processes of heat exchange that 
operate in the real atmosphere can be approximated, in a crude way, by a model 
in which the diabatic heating of the air depends only on the earth-air temperature 
difference. To test whether this simple model will capture the large-scale features 
of the heating and cooling of the real atmosphere, it is applied to the real atmosphere 
and, as far as the normal heating distribution is concerned, is found to give good 
results. The determination of the surface temperature of land having zero heat 
capacity and of ocean having infinite heat capacity is described. The importance 


of eventually adding a term that will correspond to the release of latent heat of | 


condensation when the air ascends is also discussed. 

In the fourth section the term “general circulation” is precisely defined and four 
experiments are proposed, depending on the nature of the underlying surface. In the 
first experiment the planet is entirely land-covered. In the second it is entirely ocean- 
covered, with the initial ocean temperature a function only of latitude. Because of 
zonal symmetry of the forced heating in these two experiments, they can predict 
mean zonal and mean meridional circulations but not the general circulation 
proper. Experiments of these two types have already been performed by N. A. 
Phillips, and the results of the first are summarized in this section. 

In the third experiment that is proposed the planet is again entirely ocean-covered, 
but the initial ocean temperature is a function of longitude as well as latitude. 
This experiment should predict a time-averaged motion which is not zonally sym- 
metrical and hence is a complete general circulation. An experiment of this kind is 
being programmed by Mr. Avraham Huss for the WEIZAC Computer at the 
Weizmann Institute of Science, Rehovot. Finally, an experiment is proposed in 
which the planet is partly land- and partly ocean-covered. This experiment will be 
carried out by the writer on the IBM 709 Computer at the University of California, 
Los Angeles. 

The writer wishes to thank Professors J. Charney and N. A. Phillips for the op- 


portunity to work with thei meteorology group, first at the Institute for Advanced 

Study in Princeton and later at the Massachusetts Institute of Technology, where | 
he was able to discuss and develop many of the ideas presented in this paper. He | 
also wishes to thank Mr. Avraham Huss, his assistant in the Meteorology Department - | 


of the Hebrew University, for his help in the development of these ideas. 


The writer thanks Prof. P. A. Sheppard for providing him with the data for Figure | 
Sa, and the Director of the Meteorological Office, Great Britain, for permission to 
use his manuscript maps of 500-mb normal temperatures, compiled by the Upper 
Air Climatology Branch. To Mrs. Sarah Schlesinger he is grateful for making many | 
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of the hand computations, and to Mr. Peter Grossman for preparing the illustrations 
for the printer. 


I. THE DYNAMICAL EQUATIONS FOR LARGE-SCALE FLOW IN A TWO-PARAMETER MODEL 
OF THE ATMOSPHERE 


For the study of the large-scale motion of the atmosphere we will use the vorticity 
equation in its approximate form, and in the x, y, p, t-coordinate system, 
BINVVC HN + foes ev x), a) 
where ¢ is the vertical component of the relative vorticity, V is the non-divergent 
component of the horizontal velocity, f is the Coriolis parameter and fm is the mean 
Coriolis parameter, w is individual pressure change (the “‘vertical velocity” in x, y, p- 
space), g is acceleration of gravity, ¢ is the horizontal component of the small-scale 
vertical eddy stress, p is pressure, V is the horizontal del operator on an isobaric 
surface, and k is unit vertical vector. 
The non-divergent component of the horizontal velocity and the vertical component 
of the vorticity can be written as 
v=kxVy, (2) 
¢=k-VxV=VYy, (3) 
where y is the stream function for the non-divergent horizontal velocity and V2 is 
the isobaric Laplacian operator. 


In the experiments described in this paper, we will use the quasi-geostrophic 
approximation, 


and 


& 
where Z is height and fn is the mean Coriolis parameter. 

The two parameters which will carry the history of the motion will be the stream 
functions for the non-divergent horizontal velocities at the 250- and 750-mb 
pressure levels. For convenience, we designate quantities measured at the 0-, 250-, 
500-, 750-, and 1000-mb levels by the subscripts 0, 1, 2, 3 and 4, and quantities 
measured at the surface of the earth by the subscript s, as shown in Figure 1. 


mb 
fo) 0 
! 250 
2 500 | 
3 750 
4 1000 


Figure 1 
Pressure surfaces designated by the subscripts 0, 1, 2, 3, 4 and s. 
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With this notation we then have, from (1), 


a =—YvV, ‘VC, Py) +fn(Z sie (k-Vx-) |, (5) 
d 
an oe = — V,:V(C, +f) ial?) - s[- (k-Vx7) |), (6) 


Given y; and y3 (and hence V,, V3, ¢, and ¢3) at any time ¢, and the boundary 
conditions, there are six unknowns in equations (5) and (6): dp,/0t , Cy,/ot , (Cw/dp),, 
(Aw/Op); , (Ot/ép), and (0r/ép)s « 

We replace the last four of these unknowns with more tractable parameters by 
the approximations: 


@ = @,B(p) + ee aes = Be) Vupe- (7) 


where, given w,, the first term on the right describes the distribution of w with p 
that would occur if the ground beneath were level, and the second term gives the 
deviation from w,B(p) due to the orographically forced vertical motion, V-(p/ DV s) ps ° 
(p/ps)V ups is the slope, in x,y, p-space, of an orographically determined stream surface 
which varies from Vyps at the ground to zero at the top of the atmosphere. B(p) 
is an arbitrary or empirically determined function of pressure, and Vyps is the 
horizontal gradient of the time-averaged surface pressure; 


t= C(D) Ts, (8) 
where C(p) is an arbitrary or empirically determined function of pressure, and tT; is 
the surface stress 

= Ne, (9) 
where K is the surface stress coefficient and Vs is the geostrophic surface wind; 
and 

K = Km S(v3-4), (10) 
where Km is the mean surface stress coefficient, and S(y3-4) is an arbitrary or empiri- 
cally determined function of the static stability of the lowest layer of air. 
Substituting from (7), (8), (9) and (10) into equations (5) and (6), we obtain 
as feet gim hy “VC, +f) + fin o2 mae + fin (Sv), = 


sie =), |: cape eKn(S) k- ‘VXS¥,, (11) 


R= WG +A + fan (F 2) 4 sal(2v0) 


~¥,p,(5 ae £Kn(Z) kV XSW (12) 
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Given y, and y,, (11) and (12) are now two dynamical equations in the six 
unknowns: dy,/dt, Op,/Ot, Ww, Vo, Vs and S. 
We now make use of the thermodynamic energy equation, in the form 


0 (oz Oz Rh 


where 


~ gpO- gg pe ap’ ” 
and h is the heat added to unit mass of air in unit time, T is temperature, @ is potential 
temperature, R is the gas constant, cy is specific heat at constant pressure, and 
*% = R/cypis Poisson’s constant. 
Taking the integral of (13) from Pp, to p; , with substitutions from (4) , (7) and (14), 
we obtain . 


D3 
-vs—w=-f v-v(2) Pa) =o Re 


Vups OO 2. 
B Vv B —— p*! dp. (15) 
nai “hi, () +L Vp —Y,p,Bp)] Sb = 
We now introduce the approximations: 
h = F(p) Q, (16) 


where Q(x, y, tf) is the heat added to unit column of air in unit time, and F(p) is 
an arbitrary or empirically determined function of pressure; and 


Bee =ay(p)-, for pe<pi< Py (17) 


which states that between the levels p; and p3 the lapse rate of potential temperature, 
y , is a function only of pressure, as shown in Figure 2. 


mb 
(@) 
/\ 
At y=y (x y,P,t) 
rN 
250 “ 
Pp 
: : oe a 
750 
nae yay (xy,pt) 
900 
T= 


Figure 2 
The distribution of y. 
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The approximations (16) and (17) enable us to integrate all three terms on the 
right in equation (15). 
We integrate the first term on the right in (15) by describing the wind at any level — 
between p, and p3 as 


v=V,+(V—¥)), (18) 
so that 


=f y-v(2 a he dc ‘V (v3 — —w+f (V, —V): oe ap) (19) 
Pi 
Because V e an isobaric operator, 


Cp 2 ee R R p*1 | 
SE, egret « 5 Le Ve. 20) | 
Vo - f Gp Sm ‘ Smp Sm D% ed | 
Using (17), 
Ps 
0=0;+ § ve, 
P 
or VO =VO,, forp,<p<pP3, (21) 
and substituting into (20), 
Oy sir ek 
V pte Np B: “vO, » for p, < p < D3. (22) 
For the vertical shear Ds the wind we ioe have, from (2) and (22), 
ov R p 
5p ee v= #2 B  kxVO;, for p, < Pp < D3. (23) 


In other words, the assumption that y = y(p) between p, and p, necessarily implies 
that between these two levels the wind shear, @V/dp, is constant in direction and that _ 
its magnitude is proportional to p*-1, 
Integrating (23), 
v=V, —({k x V@)) 


1 
- Pp guerre a 


=", +x V0)(F5) 7 —P9 
and from (24) and (22), 


0 R x—1 | 
(0-0 2 = — x6, V0) (ZA) = P(px—p*)=0. (25) : 


When y = 7(p) the last term on the right in (19) vanishes and 


(24) 


P3 | 
oy 
= f v (2) dp = — V;‘V(y, —,). (26) 
Pi 
Similarly, if we describe the wind at any level between p, and p; as 
vV=V, + (V—V;) (18a) 


we obtain 4 Oy 
-fovv Ci dp = —V, -V(y; —¥). (26a) 


Pi 
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From (24) we also obtain 


Pp" 
Pe fehl ae 
Ae a: 
caf +(V, ap V3) oe (27) 
P3 Pi 
and 
V2 = V; + (V; — ¥3) fete 
? Race | 1 3 Px — px 
= V3 +(¥; — V5) = eae . (28) 
Py 
Substituting now from (26), (27), (28), (17) and (16) into equation (15), we obtain 
a V 
(v1 —¥3) = —Via "Wy, — ys) + Go, +10 + (H,Y,— Hy)- 2, 29) 
Ps 
where V,,, denotes either V, or V3, ea 
G= B x1 8 30 
Fabs a (Py(p)p*"* bp, (30) 
H,=— ale (Px — p*) p*y(p) dp — 
tT Soe BE = 
— (px — p*) p, it Bp) p*“! (p) dp], (31) 
—— a. ‘ (p* — p*) p*y(p) dp — 
4 fm PS i 
— (p* — PP, f B(p) P* y(p) dp], (32) 
and 
33 
Se oe — F(p) dp. (33) 


Substituting from (27), (28) and (29) into (11) and (12) we obtain, finally, 


Ss 
a EV IG Pot) ay = Vp 2 
if Ds 


ac Gh 
— gKm (=) eV x SVs + 1,Q ’ 


and 
VuPs _ 


ABs = _V,-Y(q, +f) + (M1 — 9595)" 
t Ds (35) 


Nk — Km (F) kV XS + 1,0, | 
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where 


qo = CO = Ay, Br Y3)> (36) | 
aS C3 + A,(y, oe Y3)3 (37) | 
_ fm 
4=-2(Z). (38) 
_ fm 
= =--4 =), (39) 
=A, + are —pe(l+ )— r(2) (rps), (40) 
= A,H, — ES [ «vx + PAZ) 0% — py]; (41) 
and 
” =1,H, es [x -p(2) v9), (42) 
5 = AH, + =a li-ma +94+AZ) rp] 4) 


We note, in passing, that if we had adiabatic frictionless flow over level ground, or 
Q=0, K=0 and Vzp;=0, 
then (34) and (35) would reduce to the simple forms 


a 
Si =~ Ve (q, +i); 
and aq 


For the surface wind, V; , we write 


Ps oy 
P3 


and from (2) and (20) 


Re 
Ss eae f (k xVO) p dp. (44) 
4 P3 
Using the approximations: 
00-04-90, = 
—_S = 5 for < <, > 45 
; op Pa — Ps P3 Pp Ps ( ) 
an 


where Ty is the temperature of the ground reduced to 1000-mb; and substituting 
for VO, from the integral of (22), which is 
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75 
R P3 
Viy, — = —— WO x1 6p , 
we get é: 
vo — —*JmPa_() _ P—Ps oy _ EEE 3 48 
RE pak re e 0) + FO VT, (48) 
for p3 << p < Ps. 
Substituting from (48) into (44) , we obtain 
Vs = V3 — L(ps) (Vj — V3) — M(ps) k x VT | (49) 
- where 
Le) = ee. “lbp, (50) 
and 
M oS ae 51) 


An example of the vertical wlan : he wind from p; to ps, in this model, 
is shown in Figure 3. For a simple illustration we show only the westerly component 
of the flow, for 

um =35ms!, w=i15ms!, pp=—ps, *« = 0.286 
and the three cases, 

VT, =0 , which gives uy =11.7ms1!, 

VT, =VO; , which gives uy = 8.6ms!, 

VT, = 2V@; , which gives Us = 5.6m s-b. 


ro) 5 105 20 25 30 35 ms" 


Figure 3 
Dependence of the surface wind upon the horizontal surface temperature gradient. 
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The example shows that there can be large differences in the surface wind depending | : 
upon whether the horizontal temperature gradient at the surface is smaller, the same, | 
or larger than the horizontal temperature gradient in the layer from p; to p3. 1 

For the function B(p), which determines how the vertical velocity will vary with _ 
height over level ground, we will use | 


— p2 | 

P4 | 

That this parabolic distribution of the vertical velocity is a fair approximation | 

to the normal distribution of w in the real atmosphere is indicated by those studies | 

in which the vertical velocity is computed from the horizontal divergence of the 
observed wind. 


One such study, by Palmén (1958), is shown in Figure 4. To reduce the effect of © 
instrumental errors in the wind observations and the effect of local deviations of the 
wind due to eddies of a scale smaller than the distances between observing stations, 
Palmén averaged the divergence of the observed wind over two large areas: a western 
area (W) in which the winds in the lower levels were northerly, and an eastern area 
(E) in which the winds were southerly. These two areas and the streamlines of the 
wind at the 700-mb level are shown in Figure 4a. The vertical velocity w, computed 
from the average divergence of the observed wind in each of the two areas, is shown 
in Figure 4b. The function B =q@/w», averaged for the two curves, is shown in 
Figure 4c. 


Another study of this kind is that of Sheppard (1949), who measured the divergence 
of the observed wind over the triangular area between three radar wind stations. 
To reduce the effect of instrumental errors in the wind observations and the effect 
of eddies smaller than the distances between the observing stations, Sheppard formed 
the time-average of the divergence of the observed wind over monthly periods, 
averaging the cases separately according to the sign of the vertical velocity at 500-mb. 
The results for the triangular area (Liverpool-Larkhill-Downham Market) for the 
period 1 Dec. 1947-31 Jan. 1948 (derived from figure 3 in Sheppard’s paper and the 
number of cases used for each curve, supplied by Prof. Sheppard in a personal 
communication) are given in Figure 5a. The function B = w/w, averaged for the | 
two curves, is shown in Figure 5b. 

Both of these studies, which are the only ones of their kind known to the writer, 
indicate that (52) represents, roughly, the average vertical variation of w over level 
ground in the real atmosphere.* 

The orographically forced vertical motion, in our model, at each level is equal to | 
the scalar product of the horizontal velocity and the slope of an orographically | 


4 For the longitude-averaged vertical velocity, kinematically computed from the observed 
winds in the latitude belt 15°-30°N, cf. Vuorela (1955), figs. 9, 10. 
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Streamlines 
at 7O0Omb. 
I5 Oct. 1954 
1500 GMT 


mb 
0 

ac 
200 

~ 

~~ 
[0p 
600 
800 
1000 

see Ce Ome Sie ae S2 Get Obs i ee. 3. 3 24 © 02 04 06 08 10 
- -! 
Ww *10-4 cb s B — 
Figure 4 


The kinematically computed vertical velocity, and the function B(p), in a synoptic example 
(from Palmén, 1958). 


determined stream surface, (p/ps)Vups as shown in Figure 6. The linear variation 
of this slope, from Vyps at the ground to zero at the top of the atmosphere, was 
arbitrarily chosen. 

It should be noted that this orographically forced vertical motion differs in an 
important way from the manner in which orography is usually introduced into 
numerical forecasting models. In the usual procedure (cf. Estoque, 1957; Haltiner, 


78 YALE MINTZ Bull. Res. Counc. of Israel | 


200 


400 


—V 


600 


800 |- 


1000 : 
36: 25. - 4 <3" -3 = 0 Ie 3 8 Oo o2 04 06 O8 10 
as x107*cb s? [yas 
Figure 5 


The kinematically computed vertical velocity in the triangular area Liverpool-Larkhill-Downham 

Market during the period 1 Dec. 1947 to 31 Jan. 1948, averaged according to the sign of the vertical 

velocity at 500-mb; and the corresponding function B(p). The open circles show the average vertical 

velocity for all the cases in which the air at 500-mb was ascending (ws599 < 0), the closed circles 

the average for all the cases in which the air at 500-mb was descending (W599 > 0) (data from 
Sheppard, 1949). 


x,Y ——> 


: Figure 6 hy 
The orographically determined stream-surfaces and forced vertical velocities at the levels 1 and 3. 
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G. J. and Hesse, T. S., 1958) the orographic vertical motion at all levels depends 
only upon the slope of thie ground and the horizontal velocity at the ground. In our 
model it depends not only upon the wind at the ground but also upon the vertical 
wind shear. For example, to use an extreme case, if the wind changes from an upslope 
to a downslope direction between the lower and upper levels, the vertical velocity 
will change from ascent at the lower levels to descent at the upper levels. This oro- 
graphic component of the vertical velocity is almost exactly the same as that which 
obtains in the x, y, o-coordinate system of Phillips (1957), but it enables us to keep 


_ the other advantages of the x, y, p-system. 


For the lapse rate of potential temperature in the levels between p, and p, we will 
use the values of y given by the International Standard Atmosphere, in which, for 
Pp > 225-mb, 

y= se = A253 x pr, (53) 

Our assumption, that y=y(p) for pj) <<p<p;, is a fair approximation to 
conditions in the real atmosphere. Except in the narrow frontal zones, which are only 
about 100-mb deep, the lapse rate at these levels in the real atmosphere is in fact 
remarkably constant. 

As we have seen, the condition that y is constant in x and y signifies that between 
Pp; and p3 the wind shear, 0V/dp , does not change direction with height and that its 


_ magnitude is proportional to p*-1 . This means that the vertical wind shear, in x, y, 


_ z-space, is 
ov av 
ee ar 
4 Aas Ss eet 2) 
jad RT RO’ 
' and 
ov 1 
feos: 54 
az @”’ ) 


_ where g is the density. Because the percentage variation of @ is small between p; and p3, 
| @V/dz is nearly constant with height between p; and p3, in the model. 


| 


Near the ground, in the real atmosphere, y varies markedly in x, y and t. And in 
the model, y varies in x, y and ¢ according to (45) and (46), where both O; and T, 
are functions of x, y, t. 

In the real atmosphere y also varies markedly with x, y and ¢ in the high troposphere 
and in the stratosphere. But y for p < p, does not enter into the prediction equations 


- for our model, and therefore need not be specified. : . 
— For the function C(p), which specifies how the vertical eddy stress varies with 


height, we will assume 7 p Bx 
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or | 
w= 3p. (56) | 
ip P, | 
c— | 
0° : 
| 
| 
| 
| 
p 
| 
3 
Ps 
4 
fe) t—> t; 
Figure 7 


Vertical variations of the frictional stress, t, and of the frictional force, 0t/ Op. 


The resulting parabolic variation of the frictional stress, +, and linear variation | 
of the frictional force, 0t/dp , is of course an extremely crude approximation. 
The usual assumption made about the vertical eddy stress is that 


ov 


where pu is the eddy viscosity. 

But in our model, as we have seen, OV/0z is nearly constant with height between 
Pi and p3, so that if is taken as constant, 

ov 0 (Ov | 

Equation (57) is therefore not a suitable way of estimating the frictional stress 
when we have only a two-parameter representation of the atmosphere, and we must 
resort to some approximation such as (8). 

For the function S(y3-,), which specifies how the surface frictional force will vary 
with the stability of the air in the lowest layer, we will assume 


1.82 for y3;-, <0 
0.182 for ¥3-4 >Yp 
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where y,, is the lapse rate of potential temperature for an isothermal atmosphere and 
y = 0 is the adiabatic lapse rate. Making the surface stress coefficient vary with 
the stability in this way througha whole order of magnitude is not entirely arbitrary, 
but is guided by two observed phenomena: (1) for a given geostrophic surface wind, 
the ratio of the anemometer-measured surface wind to the geostrophic wind is about 
two to three times larger under unstable than under stable conditions (cf. Barad 
and Davidson, 1954); and (2) for a given anemometer-measured surface wind, the 
_ ratio of the surface stress to the anemometer-measured wind is about twice as large 
under unstable as under stable conditions (cf. Darbyshire and Darbyshire, 1955). 
_ Together with the usual assumption that the surface stress varies with the square 
_ of the anemometer-measured wind, these two observed conditions justify, in a rough 
_ way, the assumption that the ratio of surface stress to geostrophic wind is highly 
dependent upon the stability of the air. 


Yadiabatic Y mean Yisothermal 


; ae 


182 


0.182 
° 0.5 1.0 °C cb-! 
s. 
Figure 8 
Dependence of the surface stress on the static stability of the air. 


The dependence of the surface stress upon the stability of the air may be an 
important factor in the vorticity budget of the atmosphere. If we write 


—k-VxKVY; = —k:-VKXVs — Kk:V Xs, (59) 


| the last term always acts to decrease the absolute magnitude of the surface vorticity. 
| But the first term on the right in (59) may increase or decrease the vorticity. For example, 
_off the east coasts of the continents, in winter, a stream of cold air from the land 
| may run parallel to warmer air coming from lower latitudes, as shown in Figure 9. 
As these air currents of different temperature pass over an ocean surface which has 
a smaller horizontal gradient of temperature, there will be a gradient of Y3-4 in the 
direction normal to the wind and, according to the first term on the right in (59), 

cyclonic vorticity will be frictionally generated in the air. A simple calculation, 
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Figure 9 


Example of the variation of the surface stress with the static stability over the western oceans. 
—T is the stress of the ground upon the air, in this case producing cyclonic vorticity. 


using a moderate value of the surface geostrophic wind and a moderate gradient 
of y3-4, shows that the rate of increase of the vorticity of the air by this means will 
be of the same order of magnitude as the observed advective changes of vorticity. | 
Moreover, with the function C(p) as given, the largest part of the increase of vorticity | 
will be in the air near the ground. We will study the possible contribution of this term 
to a high frequency of east coast cyclogenesis, in the numerical experiments, by 
comparing the development of the circulation when K is allowed to vary with the 
stability, as above, with the development when K is kept constant, all other conditions 
of the experiment being kept the same. 

A word should be said justifying the use of a linear variation between surface | 
stress and surface geostrophic wind and about the neglect of the variation of surface | 
roughness: in our model. The linear relationship is used, as a first approximation, | 
because observations in the real atmosphere show that the 1atio of the anemometer- 
measured surface wind to the geostrophic wind decreases as the geostrophic speed | 
increases (cf. Gordon, 1950; Barad and Davidson, 1954). Therefore, if the stress | 
varies with the anemometer measured wind raised to the second power, it will vary | 
with the geostrophic wind raised to a power smaller than two. We are arbitrarily : 
making this power unity. 

There are also indications (cf. Barad and Davidson, 1954) that, at least over the 
smaller scales of roughness, the ratio of the anemometer-measured surface wind | 
to the geostrophic wind decreases as the surface roughness increases. This acts to 
reduce the dependence of the stress upon the roughness when the wind used is the 
geostrophic. 
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The magnitude of the average surface stress coefficient, Km, will be derived in 
Section II. 

For the function F(p), which determines how the diabatic heating is distributed 
with height, we will assume that it is constant, or 

R= Pos xX 10-2 m2t-1 , (60) 
4 

so that one-half of the total heating of unit column of air takes place in that part of the 
column which lies between p, and p,. 

Given B(p), y(p) and F(p), in (52), (53) and (60), and R = 287 kj t- dég=! c= 
1004 kj t-! deg-!, x = 0.286, fm = 10-4 s-! and g = 9.81 m s~, we obtain the 
following numerical values for the constants: 


G = 1.596 x 106 m3st-! 
H,= 5.72 x 10°m?2s-1 
H;= 11.10 x 10°m2s" 


I = 3.08 x 102m?st-1 
Ay = — 1.261 x 10-12 m-2 
A3 = 1.261 X 10-12 m-2 
By = — 253 -x10-s! 
33> > 562° K 10> st 

vy = — 3.20 x 10 sl 


z= 2.51. Xx 40-*s=! 
For the special case of level ground, or p, = p 4» we also have 


L.= 1,66 x 10-1 


' and M = 3.82 X 105 m2s-! deg. 


H 
| 


| 


\ 


The heating function, Q(x,y,t), and the expression for the surface temperature, 
T,(x,y,t), are described in Section III. 


II. EMPIRICAL DETERMINATION OF THE MEAN SURFACE STRESS COEFFICIENT, Kin 


| Equations (9) and (10) can be used to obtain the mean surface stress coefficient, 
Km, if we know the surface geostrophic wind and at the same time we compute, 
independently, the curl of the surface stress, V x ts, where the superior bar operator 


_ denotes the mean taken with respect to time. 


Vx, can be computed by using the vorticity equation (1) diagnostically. This 


_ is done by taking the time-mean and vertical integral of (1) to obtain 


Ps 
0--f V-V (64S) Op + fm Os —gk-VX%. (61) 
0 


The time interval must be chosen sufficiently long to make the mean local time rate 
of change of vorticity negligibly small. 
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Substituting for 
Ws  — & 0s Ws = — £ Os Vs ‘Vis, (62) 


where VzH; is the slope of the ground in x, y, z-space, we get 


8 
gk-Vx% =—f VeVETS) op — afm osVe-ViHh. (63) 


From the observed field of V(x, y, p, f) over a sufficiently long period of time and toa | 
great enough height in the atmosphere, we could compute, from (63), the geographical 
distribution of Vx7;, and from this and the field of ¥; obtain the geographical | 
distribution of the surface stress coefficient K{x,y) (Mintz, 1956). 

Unfortunately, the data on V are still not available in a form that will permit this 
to be done without an excessive amount of labor. But we can reduce equation (63) 
to a form which permits us, with available data, to obtain the overali mean value 
of the stress coefficient, Km . 

Using spherical coordinates, 


ad ox Othe 
k-VXt% = adore A Ty,3 — ap ees Q) : (64) 
and when the wind is non-divergent, 
PVE = aaa |" -») | 
—  @ fos ? a es ap — i cost) 
—- a (uw) + = sap mn (ce v2 ) cos b+? 30° (65) 


where A is longitude, is latitude, a is the radius of the earth, i Bi subscripts x 
and y denote respectively the eastward and northward components of the quantity. 
Also, 


ohm OsUs Os a OHs 


&fm Os Vs‘ Vas pes Se aA + Os Vs op” > (66) 
and when the surface wind is quasi-geostrophic, or Vs; = a k x(Vups), 


Sms 


gh 0s Vs 'VaHs = Fe - | 


Bea COS (2 


ar a? ae). 3 (75 a), | (67) 


Substituting (64), (65) and (67) into equation (63) and he the zonal average, } 
we obtain | 


Ce fees ) (Boos 
ag *”* p Pe eee a J uv cos2y dp — 


(— OH | 
-i pa), (68) | 
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where the superior curved bar denotes the average taken with respect to longitude 
completely around the earth, from 2 =0 to 4 = 2z. 

Integrating (68) with respect to latitude over the polar cap from @ to = 7/2, 
we obtain, for the latitude y, 


we 1 a ae LiyeOHi 
ee 2 fa 
Txs poicns le aap 4 uv cos? p dp eae (7 aa) (69) 
In general, 
a =u = uv +uy, (70) 


where the prime superscript denotes the instantaneous local deviation from the 
longitude-mean. But because V is the non-divergent component of the horizontal 
velocity, (VY=k x Vy), we have 


3 
” Dna cos p s aa oe oe 
and 0 
uw =u. (72) 


Substituting from (72) into (69), multiplying through by 27a3 cos? » , and integ- 
rating once again with respect to latitude over the polar cap, from g to p =a/2, 
we obtain 


m2 Ps aw 
= 22a? TF 
2na3 f T2,s COS*p dg = ore" cos? y f uv’ dp— 
? & 0 
eH 
= s 
2D? f (5 7 ) cosy og. (73a) 


? 
For convenience of discussion, we will denote the three terms in equation (73a) 
respectively as Ty , Jo and Pa», or 


Tz =Jq — Po. (73b) 
| To is the net torque of the atmosphere upon the earth due to surface friction in the 
| polar cap bounded by latitude gp. 
| Pg is the net torque of the atmosphere upon the earth in the polar cap bounded by 
|| latitude g due to zonal pressure differences across the orographic features of the 
' earth’s surface; it is positive when there is a positive correlation between surface 
air pressure and the zonal slope of the earth’s surface, as illustrated schematically in 
Figure 10b. 

_ Jog is the horizontal eddy flux of angular momentum into the polar cap through the 
Jatitude surface m. It is positive if there is a positive correlation between the wind 
- components u’, v’, as illustrated schematically (for a westerly wave) in Figure 10c. 

If, instead of the approximate vorticity equation (1), we had begun with the complete 
vorticity equation and integrated over the polar cap, as above—or if we had begun 
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Pole 


Figure 10 


a) Angular momentum budget of the air in the polar cap bounded by latitude g = 35°. 
b) Distribution of pressure which produces a westerly torque upon the earth. 
c) A westerly wave which transports momentum polewards. 


more directly with the primitive equation of zonal motion—the final expression for 
tT: corresponding to (73), would contain some additional small terms (cf. Mintz, 
1955a). The largest of these additional small terms is 22a? cos? (” a Opie: 
which is the poleward flux of angular momentum through the latitude surface » due 
to the mean meridional circulation F. But latitude 35° is near the position which 
separates the positive Hadley cell from the negative Ferrel cell of the mean meridional 
circulation, and the momentum flux due to the meridional circulation is therefore 
very small at this latitude and can be neglected (cf. Mintz and Lang, 1955; Bjerknes 
and Venkateswaran, 1957). 

The true role of the poleward horizontal eddy flux of angular momentum in the 
maintenance of the mean zonal motion of the atmosphere was first recognized by 


| 


i 
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Jeffreys (1926, 1927, 1933). He showed that the previously held view — that the 
large-scale atmospheric disturbances were perturbations upon a basic zonal current 
which was capable of an independent existence — was incorrect, and that the mean 
zonal wind could not exist without the eddies. But it was not until many years later 
that the order of magnitude of the eddy flux of angular momentum was first evaluated 
from the aerological data: by Priestly (1949), using the radar-observed wind to high 
levels at a single station in southern England; by Widger (1949), independently, 
using the geostrophic wind at sea level, 700- and 500-mb around the northern 
hemisphere; and by Mintz (1949), independently, using the radar-observed wind to 
high levels at 10 stations across North America and the 500-mb geostrophic wind 
around the northern hemisphere. 

In addition to the magnitude of the momentum flux, these and later studies have 
shown that the ratio of the large-scale eddy flux of angular momentum to the space- 
averaged or time-averaged gradient of angular momentum is not constant, but varies 
in space and time. In fact, the flux is sometimes observed to be in the direction opposite 
to the gradient (so that the fluid may then facetiously be described as having a “‘nega- 
tive viscosity’’). Moreover, the observed ratio has no apparent relation to the ratio 
of flux to gradient of other atmospheric properties, such as the sensible heat. For 
these reasons the “Grossaustauch”’ concept of Defant (1921), which by analogy with 
molecular and small-scale eddy motions describes the maintenance of the mean 
zonal wind as a large-scale diffusion process, .is incorrect. We are dealing with a 
thermally active system in which over one part of the spectrum of motions there is a 
transfer of kinetic energy from the smaller scale to the larger scale motion. Therefore, 
to understand the dynamics of this system, we must integrate the differential equations 
which describe the individual details of the velocity in that part of the spectrum of 
motions which generates the kinetic energy. It has been known for some time (i.e. 
Bjerknes, 1935; Mintz, 1947) that the kinetic energy of the general circulation is 
generated by the cyclonic-scale motions—by the familiar cyclones and anti-cyclones 
of the weather map and their associated upper level circulations, which have a 
characteristic wave length of 3000 to 4000 km. 


To insure that the data for all of the terms in our computation of the surface 
stress coefficient shall be homogeneous, we will use the numerical data from the 
U.C.L.A. general circulation investigation. 


wn 


Figure 11 shows the quantity u’v’ at latitude 35°N, for the months January- 
February-July-August 1949, from the U.C.L.A.study (Mintz, 1955b). The winds used 
are geostrophic, with the approximations 


| 1 AY 1 AY 
gates cats See eae. 
4 a Ay Bea). COS pm fA 
where : aie ne 
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u’y’, computed from the geostrophic winds at 35°N, in Jan.-Feb.-July-Aug. 1949. 


Using these values of uy’ , and integrating from 1000- to 100-mb—thereby 
neglecting the flux in the uppermost ten per cent of the mass of the atmosphere, 
which extrapolation of the upper part of the curve in the figure suggests would be 
only a small error—we obtain 


J35° = 2.86 X 10!5t m2 5-2. (75) 


It is important to compare the momentum flux obtained from the geostrophic wind 
in this way with the flux computed from radar wind and pilot balloon observations. 
When the geostrophic wind is computed by taking contour height differences over ten 
degrees of latitude and longitude, as above, which at 35° latitude are distances of 
about 1000 km, we are, in effect, smoothing out all disturbances in the wind field 
of a smaller size and are computing the momentum flux due only to eddies having a 
wave length longer than 1000 km. On the other hand, the wind as measured by the 
standard pilot balloon or radar tracking method is the wind averaged over a two 
minute interval during which the balloon drifts horizontally only some few kilometers. 
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The momentum flux computed from these winds will therefore contain the contribution 
of all horizontal eddies larger than a few kilometers in wave length. 

Fortunately, at latitude 30°N we have momentum flux computations ona sea- 
sonal basis computed from both kinds of winds. The flux at 30°N computed 
from the geostrophic wind for the two summer months, July-August 1949, was 
1.27x10!© tm?s-2, and for the two winter months, January-February 1949, it 
was 4.50 x 10!© t m?s~? (Mintz, 1955b). In addition, we have the poleward eddy 
flux of angular momentum computed from the radar and pilot balloon observed 
wind at 19 stations distributed around the globe in the vicinity of latitude 30°N 
(Starr and White, 1952): for the six summer months, May to October 1949, this flux 
averaged 2.05 10!©t m2s-2, and for the six months. February, March, April, 
November and December 1949 together with January 1950, considered as forming a 
composite winter season, it averaged 3.95 x 101° t m2 s-2. 

Assuming that the flux varies with only an annual period and with its maximum in 
mid-winter, these values give us, respectively, 

Japs peaate: = (2.9 + 1.7sin©) x 10!6t m2s-2, (76) 


ae J39°, rawin-pibal = (3.0 + 1.5 sin ©) x 10!6t m2 s-2, (77) 


where © denotes the time of the year (© = 0, 2 on November 1, © =2/2 on Febru- 
ary 1, etc.) (cf. Figure 12). 


5 
x10'° tm?s-? 


Nov Dec Jan Feb Mar Apr May Jun July Aug Spt Oct 


Figure je 
The poleward eddy flux of angular momentum at latitude 30°N, computed from the radar and pilot 
balloon winds (Starr and White, 1952) and from the geostrophic winds (Mintz, 1952). 


The fact that the two fluxes are so very close in mean magnitude and in seasonal 
amplitude is probably to some extent coincidental. What is significant is that the 
flux computed from the radar and pilot balloon wind is not substantially larger 
than the flux computed from the space smoothed geostrophic wind. (The same is 
also true in other latitudes, where we can compare the mean annual fluxes, and where 
there is a large gradient of the mean zonal wind: cf. Starr and White, 1955; Mintz, 
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1955b). The interpretation we must put on this is that eddies in the range of about | 
1 to 1000 km in wave-length contribute very little, if anything, to the horizontal 
flux of momentum. If we assume that eddies smaller than 1 km also contribute very | 
little to the horizontal momentum flux (cf. Panofsky and Crutcher, 1951), then | 
the lateral eddy stress can be neglected compared to the advective transport of | 
momentum by the cyclonic-scale geostrophic wind. The agreement between (76) | 
and (77) justifies the approximation made in equation (1), where the frictional force 
is dependent only on the vertical eddy stress and the lateral eddy stress is neglected. 


Before we can compute P, , the mean westerly torque upon the earth due to the 
correlation between surface air pressure and the slope of the earth’s surface, we must 
first decide what part of the spectrum of the earth’s topographic irregularities we 
will include in this term. It is clear that the dividing line between this term and the 
Tz» term is an arbitrary one, because what we call the surface (or skin-friction) stress 
may include the forces due to horizontal pressure differences across individual hills, 
individual trees and blades of grass, or the horizontal pressure differences across 
individual grains of sand. 

The distinction between 2. and Té must be an operational one, depending upon the 
way in which we intend to use the prediction equations (34) and (35). In the experi- 
ments described in this paper these equations will be replaced by their finite difference 
equivalents in which all gradients, including Vyps, the slope of the earth’s surface in 
x,y,z-space, will be replaced by differences over two grid intervals. The grid interval 
we will use will be of the order of 200-400 km. Therefore, the surface stress coefficient 
which we require is one which will include the horizoatal forces due to all scales of 
surface roughness up to twice this grid interval, leaving for the orographic term in 
our computations only that part of the spectrum of the earth’s irregularities which 
is larger than two grid intervals. 


Accordingly, for the computation of the Po term we have used the space-averaged - 
topography given by Berkofsky (1955), in which the averaging was done over intervals 
of 5° of latitude and longitude. 


For convenience of computation, we let 


P=ptpD 
Then = oH _ 20H, aH 
and 
OH, ae Peete, 
> a) = G ay’ v8) 


where the prime superscript denotes the deviation from the longitude mean. 


‘By linear interpolation between p’ at sea level and p’ at 10,000 feet (obtained from 
the 700-mb height field), averaged for the four months Jan.-Feb.-July-Aug. 1949, 
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p was obtained for the ground level given by the space smoothed topography. From 
this and 0H;/0A given by the space smoothed topography, we obtain the distribution of 


[7% OH,/0A) 1g shown in Figure 13. 


35°N 


2.0 xt0*t s-? 


1.6 


Figure 13 
[= pea rte Ne 


The pressure term, Lp’ (0 A:/0 a). , for Jan.—Feb.—July—Aug. 1949. 


Integrating from latitude 35°N to the pole, we obtain 
P35° = 0.67 x 10!5t m2s-2. (79) 


White (1949) computed this term from maps of normal pressure, using a topography 
- broken down into 1-km vertical blocks, and (neglecting the torque north of 65°) 
_ obtained, for the yearly mean, P35° = 0.59 x 10!6t m2s-2. Of course, if we were 
to use a more detailed surface topography we would obtain a larger value for this 
term and consequently a smaller value for the surface stress coefficient. But this is 
- precisely what we do not want for the experiments we are now concerned with. 

_ From (73), (75) and (79), the mean zonal torque of the atmosphere upon the earth, 
| due to surface friction in the polar cap north of latitude 35°, is therefore 


T35° = (2.86 — 0.67) x 10!6t m2s-? = 2.19 x 10!6t m2s-2. (80) 


_ From equations (9) and (10), 
fi: 7/2 7/2 


: . 2a f Tins cos2p ég =2na3 Kn f usS(73-4) cos? dp. 
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Neglecting any correlation between the geostrophic surface wind and the static | 


stability, this gives us 
To 


i a (81) | 
UgS(3-4)m 
where " 2 | 
Up = 2na3 { Us COS*p Op . (82) | 
? 


The latitudinal distribution of the mean zonal geostrophic surface wind, u;, for 


the four months of 1949 is shown in Figure 14. The zonal surface wind, u; , was 
first found for the surface elevation given by the space smoothed topography by 
using linear interpolation between the wind at sea level and the wind at 700-mb 


(10,000 feet), at every 5° of longitude. Then the longitude average of u; was taken to 
obtain us; . 


35°N 


Figure 14 
The geostrophic zonal wind at the earth’s surface, sea level and 700-mb, averaged over all longitudes 
; in Jan.—Feb._July—Aug. 1949. 


For the purpose of comparison, Figure 14 also shows as and as , the zonal 
winds averaged at sea level and at 700-mb individually. We note that at sea level, in this 
period, the transition from tropical easterly to mid-latitude westerly zonal wind was 
at 43°N, and at 700-mb the transition was at 19°N. But the mean surface zonal wind 


| 


( 


changed from easterly to westerly at 32°N, which is almost exactly the latitude of the | 
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observed maximum in the poleward angular momentum flux for this period (Mintz 
1955b). 


Using the observed distribution of iley we obtain, by integration of (82), 


U35° = 6.90 x 1020 m4 5-1, (83) 
The mean value of y,_, is 0.5°C cb-!, and from (58), 
S(y3-4)p, = 1.0. (84) 
Substituting, finally, from (80), (83) and (84) into (81), we obtain 
Km = 0.32 x 10-4 t m-?s-! = 0.32 dyne cm-2/m s-1. (85) 


This value for the surface stress as a function of the surface geostrophic wind is 
_ about 14 times as large as that which can be derived from the conventional ways of 
estimating the surface stress. Phillips (1956), for example, obtained and used the 
_ equivalent of Km = 0.20 x 10-4t m~? s-!. 
_ But if we were to neglect the pressure term in (73), so that every scale of surface 
| roughness was included as part of the skin-friction stress, then, from (80) and (81), 
| we would obtain Km = 0.42 x 10-4t m2 s-!, 


| Ill. THE HEATING FUNCTION, Q(x, y,t) 


We can characterize the atmosphere-earth system as a thin layer of gas gravitationally 
| held to a rotating sphere upon which there falls a beam of parallel radiation from 
a distant source, the declination of the beam varying with an annual period. Very 
little of the incident radiation is absorbed by the gas, but a sizable part is reflected 
and the remainder is absorbed by the underlying sphere. By involved and complex 
_ processes of heat exchange there is a net transfer of energy from the sphere to the gas. 
_ At the same time, the gas itself radiates energy to space. 
The net heating of a unit column of the gas, or atmosphere, Q, is equal to the heat 
_ transfer from the earth to the atmosphere, Q4, minus the heat radiated by the atmos- 
_ phere to space, Qo. 
For the energy radiated by the atmosphere to space, Qo, we will assume, in our 
- model, that all of it originates at some level, p < p; , at which the horizontal tempe- 
_rature gradient is negligibly small. (Up to this point, we have not specified in any 
_ way the temperature distribution above the level p;). Consequently we assume, in the 
| model, that the energy radiated to space, Qo, is a constant. This treatment of the 
outgoing atmospheric radiation is very similar to that used by Simpson (1929) for 
the real atmosphere. 
In the real atmosphere, the heat transfer from the earth to the air takes place in 
-a complex way. Part of the transfer is due to radiative exchange between earth and 
| air, part is due to conduction at the earth-air interface, and part is due to evaporation 
of water from the earth (cooling the earth) and condensation in the atmosphere 


(warming the air). 
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In the real atmosphere, the heat transfer by radiation depends on the ground-air | 


temperature contrast, the water vapor distribution and the cloudiness. The heat 


transfer by conduction depends directly on the small temperature difference at the | 
earth-air interface, but indirectly it depends on the temperatures above the interface | 
because the temperature lapse rate in the lower levels of the air determines how | 
rapidly eddy diffusion transfers the conducted heat away from the interface. And the — 
heat transfer by evaporation-condensation depends to a large extent on the tempera- | 
ture lapse rate in the lower layers of the air because much of the condensation and | 


precipitation is of the convective type, especially in the lower latitudes. 


In our model we will make the heat transfer from earth to atmosphere, Q4, — 
proportional to the temperature difference between the ground and the air —the 


temperature of the air, for this purpose, being taken as the temperature at the level 
Pz. Because we have set Ty ~ T4, this is equivalent to making the heat transfer 


from earth to air dependent on the temperature lapse rate in the lower half of the | 


atmosphere. 
With these assumptions, we have 


O(x,y,t) = Q4(%,y,t) — Qo = DAT (x,y,t) — Qo (86) 
where AT = (7, —T>) and b is the constant of proportionality. 
Taking the global average of (86), we obtain 
Om = D(AT)m — Qo = 9, (87) 


where the subscript m denotes the time- and space-average over the entire system. 
Substituting for Qo from (87) into (86), we obtain the simple heating function 


OQ = bIAT — (AT) mn]. (88a) 


According to this heating function, whenever the lapse rate in the lower half of 
the atmosphere is larger than the mean lapse rate the heat transfer from earth to 
air will be larger than from the air to space, and there will be a net heating of the 
atmospheric column. But when AT is smaller than (AT)m, the upward heat flow 
from earth to air will be smaller than the heat loss to space and there will be a net 
cooling of the column. These two conditions are illustrated schematically in Figure 15. 

For (AT)m we will use the mean 500-mb temperature, T),m = — 15.2°C (Burdecki, 
1955), and the mean surface temperature, Ty,m = 15.6°C, obtained from the normal 
sea surface temperatures (Sverdrup, 1942) and, over the lands, the normal surface 
air temperatures reduced to sea level (Hann-Suring, 1939). These give us 


(AT) m = 30.8°C. (89) 


| 


= 
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Figure 15 
The net heating of an atmospheric co umn as a function of the earth-air temperature difference. 
AT = Ty —7>. 

For Qo we will use the average value, Qo,m, from recent estimates made for the 
real atmosphere by different authors. In the real atmosphere, unlike our model, 
there is some direct absorption of solar radiation by the atmosphere. To obtain 
the net radiative exchange that takes place directly between the atmosphere and 
space, we will subtract this direct absorption from the long wave radiation emitted 
by the atmosphere to space. The estimates are shown in Table I. 


TABLE I 
Mean radiation, in 10° cal m~? day 


Alissow, Drosdow, 


London Houghton and Rubinstein 
(1951) (1954) (1956) 
Long wave radiation emitted by theatmospheretospace 41 41 38 
Short wave radiation absorbed by the atmosphere 9 13 14 
Difference 32 28 24 


Taking the average of these values, we have 


Qo = 28 x 105 cal m™ day!. (90) 
Substituting from (89) and (90) into (87), we obtain 
b =0.90 x 10° cal m-? day! deg™}, (91) 


_and from (88a), 


QO =0.90 (AT — 30.8) x 10° cal m™? day-'. 
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As a check on whether this heating function is able to reproduce the essential | 
features of the heating and cooling of the real atmosphere, we have applied (88b) | 


to the real atmosphere. For the temperature, T,, we used the world charts of normal 
500-mb temperature, for the winter and summer seasons respectively, compiled 
by Upper Air Climatology Branch of the Meteorological Office of Great Britain. 
These charts were kindly supplied to the writer in manuscript form by the Director 


of the Meteorological Office. The winter and summer surface temperatures were | 


taken, as before, from Sverdrup and Hann-Suring. The results of the test are shown 
in Figures 16 to 18. 


NP 70°60° 50° 40° 30° 20° 1o° Eq. 0° 20° 30° 40° 50° 60° 70°S.P 


x10°cal m-*day-! 


Observed eddy heat flux 


x!0"cal day-' 


41, ' 
NP 70° 60° 50° 40° 30° 20° 10° Eq. 10° 20° 30° 40° 50° 60°70° oe 
bd 
: Figure 18 
(a, b, c) The computed heating of the atmosphere in January, July, and for the year, as a function 
of latitude. 


(d) The atmospheric poleward heat flux. The continuous line shows the computed total fiux 
The broken line shows the observed eddy flux. 
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During northern hemisphere winter (Figure 16), our heating function shows a net 
cooling of the atmosphere over Canada, the United States and Greenland and over 
central and eastern Europe and Asia, with the maximum cooling in northeastern 
Asia. A broad band of heating extends from southwest to northeast over the At- 
lantic, with its maximum between Norway and Greenland. There is a heating max- 
imum over the Mediterranean and another maximum off the east coast of Japan. 
A band of heating extends from the central south Pacific, across Australia, where it 
has its maximum, and across the south Indian Ocean to Africa and to the eastern 
Atlantic north of the equator. Cooling is found in the eastern parts of the subtropical 
southern hemisphere oceans. In the middle and higher latitudes of the southern 
hemisphere there is cooling. 

During the northern hemisphere summer (Figure 17), a band of heating extends 
across the tropical Pacific, across central Asia, with a maximum over the highlands 
of Asia, across southern Europe and north Africa, with a relative minimum over the 
eastern Mediterranean and a strong maximum over Africa, across the tropical Atlantic 
and across western North America, where there is another maximum over the 
land of high elevation. There is cooling in the higher latitudes over the continents, 
and cooling in middle latitudes over the oceans and over the eastern parts of the 
subtropical oceans. (Cooling signifies a smaller upward heat transfer from earth to 
air than from air to space.) In the southern hemisphere, in its winter season, there 
is a middle latitude band of heating circling the globe corresponding to the band 
of heating over the north Atlantic during northern hemisphere winter). Between 
the southern hemisphere middle latitude band of heating and the northern 
hemisphere tropical band of heating, there is a zone of cooling, with maxima 
over Australia, and the eastern south Atlantic and Pacific. In the higher southern 
latitudes, there is cooling. 

Of course, all that we have shown in these two figures is the observed distribution 
of AT, the vertical temperature gradient in the lower half of the atmosphere, express- 
- ed as the net heating of the atmosphere according to (88b). But it is surprising how 
- well these patterns conform to our preconception of what the true distribution of the 
| heating and its changes with the seasons ought to be. Unfortunately, we do not 
_ know the quantitative geographical distribution of the true heating, so that we cannot 
make any more critical comment on Figures 16 and 17 than to say that the 
| patterns “look reasonable”’. 

A more critical test of the heating function is obtained by taking the zonal average 
of the heating. The zonal averages of the fields given in Figures 16 and 17 are shown, 
respectively, in Figures 18a and 18b (poleward of 60°S, the temperatures T> 
- were taken from Burdecki, 1955). The figures show net heating in the lower latitudes, 
net cooling in the higher latitudes. In January, there is zero net heating at 37N and 
' averaged over the whole of the northern hemisphere there is a net cooling (which 
- can be balanced by a northward heat flux across the equator). In July, the latitude 
of zero heating has shifted to 58°N and averaged over the whole of the northern 
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hemisphere there is a net heating (which can be balanced by a southward heat flux 
across the equator). It will be noticed that averaged from pole to pole there is a very 
small net cooling in January and net heating in July. 

One shortcoming in the distribution of heating shown in these figures is the absence 
(or underestimate) of the narrow belt of intense heating which we expect to find in the 
regions of the heavy intertropical convergence rains. But this, it should be pointed out, 
is very likely not a shortcoming of our heating function but of the temperature data 
used in the computations. The 500-mb temperatures that we used showed broad 
maxima of temperature in the equatorial regions. But the detailed temperature analysis 
over the Atlantic by Schnapauff (1936) has shown that the air in the middle troposphere 
is colder along the line of intertropical convergence than it is on either side. This 
is probably true for the line of intertropical convergence generally, and if we were to 
introduce such a minimum in the 500-mb equatorial temperature field our heating 
function would show the expected equatorial heating maxima. 


Figure 18c shows the average of the January and July computed heating. Integrating 
this heating over the polar caps bounded by each latitude , we obtain the poleward 
atmospheric heat flux shown by the continuous line in Figure 18d. Beginning at 
about latitude 4°S there is a heat flux directed toward each pole, with the two maxima, 
4.7 x 10!9 cal day! at 43°N and 4.3 x 1019 cal day—! at 51°S. 


For comparison we show the observed northern hemisphere eddy flux of heat 
(Mintz, 1955b) by the broken line in Figure 18d. As we pointed out earlier, 35°N is 
about the latitude which separates the Hadley cell from the Ferrel cell of the mean 
meridional circulation. Therefore the heat flux due to the mean meridional circulation 
(like the momentum flux) is very small near this latitude. We notice, in Figure 18d, 
that the curves for the computed total flux and for the observed eddy flux cross at 
latitude 38°N. In the lower latitudes, where the Hadley cell operates, we expect the 
total flux to be greater than the eddy flux, and this is what the figure shows. North 
of 35°, where the Ferrel cell operates, we expect the total flux to be smaller than the 
eddy flux, and again, with just a small discrepancy, our computations show this. 
Near the pole, north of about 70°, there may be a second direct mean meridional 
circulation: but this second direct cell, if it exists at all, is very weak and will 
produce only a negligible heat flux. 


The remarkably accurate poleward heat flux given by our heating function, as 
shown in Figure 18d, must in part be considered coincidental: so good a result 
obtained from so simple a model must in part be due to a fortuitous cancellation of 
errors. But though we may doubt that the heating function will be this good in all 
respects, the results of the test should encourage us to believe that we have a useable 
first approximation. 


When the heating function (88), is used for numerical prediction experiments, 
it is necessary also that the surface temperature, T,, be determined. This we will do 
for the ocean and land areas separately in the following simple way. 


ee) 
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For the oceans we will assume in the model an upper layer of water in which 
there is no vertical temperature gradient and which has no heat exchange with the 
water below it. For this layer the local change of temperature can be written 


OT, Ss =0; 
Bypse V-VT7+ G (92) 


where S is the solar radiation absorbed by the ocean and C;, is the heat capacity of 
unit column of the upper layer. 

We willassume also that in the ocean the advective temperature change is negligibly 
small and that the heat capacity of unit column is very large, so that 


and 
| Ty =(Ty)t=o, Over the ocean. | CY 


Over the land we will assume in the model that the temperature changes of the 
ground are confined to a very shallow layer having a heat capacity which is negligibly 
small. Then 


Ch ot = S— Q4~ 0, (94) 
and 
S='0, (95) 
or LS LE LO AOE LE 
S 
Ty = oe + T> » Over the land. (96) 


The constant ocean surface temperature given by (93) approximates the observed 
condition in the real oceans in which, averaged over 5° latitude-longitude intervals, 
the local temperature changes are less than a degree centigrade per month. 

Over the land, our model is not as good an approximation to real conditions as 
it is over the oceans. But it does reproduce one of the important features of the real 


- atmosphere, which is that the temperature difference (T, — T) varies latitudinally 


and seasonally in the same sense as the insolation, S. In low latitudes and in the 


“summer season the temperature lapse rate in the lower part of the atmosphere is 


large; in high latitudes and in the winter season it is small. This is observed in the 
real atmosphere only over the lands and not over the oceans: a distinction which is 


preserved in our model. 
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It is possible to construct a more refined heating function than the one described, 
even for a two-parameter model of the atmosphere (cf. Mintz, 1958). One of the 
modifications that we may experiment with consists of treating the heat of conden- 
sation in two parts: a convective condensation which is made a function of the 
ground-air temperature difference and is part of the heat transfer, Q4, as before; 
and a non-convective condensation which is made a function of the large-scale 
vertical velocity, w, by assuming that some predetermined fraction of the air column 
is saturated wherever the air ascends. (We would not wish, at this stage, to add another 
parameter to the model in order to carry the history of the water vapor content of 
the air explicitly.) 

We can anticipate a number of important consequences of having the air in the 
model release latent heat of condensation when it ascends. One result would be the 
same as if there were an overall reduction in the stability, y; namely, a decrease in 
the wave-length of the disturbances of maximum growth rate. But unlike the effect 
of an overall reduction in y, having the air descend dry-adiabatically and ascend 
with individual temperature changes less than the dry adiabatic will produce an 
asymmetry in the disturbances that develop, smaller areas of more intense cyclonic 
development being balanced by larger areas of less intense anti-cyclonic development. 
But the most important effect of all would be that for the wave-lengths of maximum 
growth rate a smaller vertical shear would be required for a given rate of growth, 
and therefore smaller vertical shears (and smaller horizontal temperature gradients) 
would develop in the model. 


IV. FOUR NUMERICAL EXPERIMENTS 


We define the general circulation as the time-averaged velocity field of the atmosphere 
as seen in its world-wide distribution*. By the velocity we mean here the vertical 
velocity, w , as well as the horizontal velocity, V , where, as in the previous chapters, 
these are not physical-point velocities but velocities averaged over horizontal areas 
of the order of 500 km in diameter. 

For the specification of how the time average is to be taken, we define both a 
“transient general circulation’ and a “normal general circulation’. 

By the transient general circulation we mean the velocity field averaged over a time 
interval which is an order of magnitude larger than the characteristic period of the 
synoptic-scale atmospheric disturbances. As these disturbances, the ordinary cyclones 
and anti-cyclones of the weather map, have an average period of about three days, 
the transient general circulation will be the velocity field averaged for about a month. 
Thus, the transient general circulation for 15 January 1959 would be the world-wide 
fields of V and w averaged from 1 January to 31 January 1959. 


* Some meteorologists use the term “general circulation” for what we call the ‘‘mean zonal cir- 
culation”, or for the “mean meridional circulation” (or for both). Others give still different meanings 
to the expression “‘general circulation’’. It is to avoid ambiguity that we give here our own, perh : 
cumbersome but necessary, definitions. eaten 


2 
4 


Vol. 7G, 1958 YALE MINTZ 103 


By the normal general circulation we mean the transient general circulation of 
a particular day of the year averaged over a large number of years. Thus, the normal 
general circulation for 15 January would be the world-wide fields of Vand w averaged 
from 1 January to 31 January for many years. In this way we can speak of a normal 
general circulation for each day of the year. 

In addition to the general circulation proper, we can take the zonal average of the 
general circulation (the average over all longitudes). The zonal average of the west- 


east (or zonal) component of the time-averaged velocity, #, is called the mean zonal 
circulation. The zonal average of either the south-north (or horizontal meridional) 


component of the time-averaged velocity, V. 7, taken alone, or both V 7 and @ Fr together 
(which are related in a simple way by the equation of continuity of mass) is called 
the mean meridional circulation. The mean zonal circulation and the mean meridional 
circulation are functions of only two of the space coordinates, latitude and elevation, 
and are usually shown in meridional cross sections. 

Finally, we note that we can have a transient or a normal mean zonal circulation, 
also a transient or a normal mean meridional circulation, with the time-averaging 
in each case being specified in the same way as for the general circulation proper. 

In the pages that follow we will propose four numerical experiments which deal 
with either the normal general circulation or with the normal mean zonal and mean 
meridional circulations. In the first experiment the prediction equations (34) and (35) 
and the heating function (88) are applied to a planet which is entirely land-covered. 
In the second experiment these equations are applied to a planet which is entirely 
ocean-covered and in which the initial surface temperature is a function only of latitude. 
In the third experiment the equations are applied to a planet which is entirely ocean- 
covered and in which the initial surface temperature is a function of both latitude and 
longitude. In the fourth experiment the equations are applied to a planet which is 
partly covered by land and partly by ocean. 


1. The general circulation on a land-covered planet 


On a land-covered planet we have, from (95), 


Q4=S8 


and from (86), 


ORS eels OF (97) 


If we neglect cloudiness and variations of surface albedo, then Sis a known function 


of latitude, y, and the time of the year, ©. And if, for the the purpose of a simple 
_ experiment, we neglect also the seasonal variation of S, then the net heating per unit 


atmospheric column, Q, is a function only of latitude, or 


Q = S(v) — Qo = OY). (98) 
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Figure 19b shows the mean annual distribution of the net heating, O(~), as obtained 
by making S(y) proportional to the mean annual insolation and with its global 
average Sm = Qy = 28x 105 cal m~? day. This, of course, is not exactly the same 
as the heating computed from the earth-air temperature difference shown in Figure 
18c. But they are approximately the same. 

Heating which is a function only of latitude was used by N. Phillips (1956) in his 
remarkable pioneering experiment in the prediction of the normal mean zonal and 
mean meridional circulations of the atmosphere. Although many meteorologists 
have contributed greatly to our understanding of the physical processes by which 
the mean zonal and meridional circulations are maintained, Phillips’ experiment was 
the first successful attempt to predict these circulations. It was the first time that the 
observed circulations were predicted by the direct application of the laws of fluid 
motion, starting with an atmosphere at rest. 

Because we can interpret Phillips’ experiment as being the first in the series of four 
experiments which differ from one another only in the kind of surface on which the 


190 


15 x10° 
cal m-?day-' 


Q , Phillips (1956) 


10,000 Km 5,000 () 
\¢ 
5 xlo5 
cal m-2 day-! 
10 
5 
0 Q 
=) 
-10 
‘ ‘ 15 
Pole 60° 50° 40° 30° 20° (Oe Eq. 
P 
Figure 19 


a) The heating function used by Phillips (1956). 
b) The mean annual heating of the real atmosphere. 
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atmosphere rests, we will briefly describe some of its features. For details the reader 
is referred to Phillips’ paper. 

Phillips considered the motion that would develop in a two-layer model of the 
atmosphere in a rectangular region bounded by “walls” at the north and south, but 
open at the western and eastern sides (x = 0 and L), with the condition imposed 
that the flow at x = 0 is identical with that at x = L (cyclic continuity). The distance 
between the north and south walls was taken as 10,000 km (approximately the pole- 
to-equator distance on the earth). L was taken as 6,000 km. The heating was made a 
linear function of the south-north coordinate, y. The heating used, in units of cal. 
per meter square column per day, is shown in Figure 19a. 

The dynamical prediction equations were derived by Phillips for a discontinuous 
two-layer atmosphere rather than for a two-parameter representation of a continuous 
atmospehre, as in our Section I. However, except for the choice of constants, the 
resulting prediction equations are the same 

Phillips assumed flow over level ground so that, in terms of our equations (34) 
and (35), Ps =p4 and V,,ps; = 0. For C(p) he used the equivalent of C =(p/p) — 1 
for po < p < ps4, and C =0 for p < py, so that surface friction operated on V3 but 
had no direct affect on V;. For the equivalent of our Km he used 0.2 x 1074 t m=? s7}, 
and S(y3-4) = S(y3-4)m = 1. For V;= V4 he assumed a linear extrapolation of the 
wind from VY, and V3, so that in terms of our equation (49), 7 =0.5 and M=0. 
The boundary conditions used were Vy = 0 at the north and south walls. For 7 he 
chose 43 = — a, = 1.5 X 10712 m™”. In addition, Phillips used a small-scale 
lateral eddy diffusion term both for momentum and for heat, with the lateral diffu- 
sion coefficients A, = A; = 10° m’s~!. Butitis not likely that these played anything 
but a minor role in determining the form of the motion which developed. 

Figures 20, 21 and 22 show the predicted fields of the mean zonal and mean me- 
ridional circulations and of the mean temperature (for the 10th to 20th day of the pe- 
riod when disturbances were allowed to develop), together with the mean annual 
observed fields of these parameters for the real atmosphere. 

For the zonal circulation (Figure 20), the prediction given by Phillips’ experiment 
shows in a remarkable way the essential features of the real atmosphere; that is, 
surface easterlies in the low and high latitudes and surface westerlies in middle la- 
titudes, with a general increase of westerlies with elevation. (The model, in our in- 
terpretation, makes no prediction of the temperature above the level p; and there- 
fore no prediction of the wind above this level.) 

It is clear from the prediction equations that where there are surface easterlies 
there must be a net divergence in the poleward eddy flux of momentum and where 
there are surface westerlies there must be a net convergence in the poleward eddy 
~ flux of momentum. The numerical values for these flux divergences, at levels 1 and 3, 
are given by Phillips in his paper. 

One noticeable difference between the predicted mean zonal circulation and the 
observed circulation, shown in Figure 20b (Mintz, 1954), is that the predicted ver- 
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a) The predicted mean zonal circulation, in m sl (Phillips, 1956). 
b) The observed mean zonal circulation, inms ~. 


tical shear of the zonal wind, and hence the magnitude of the high level zonal wind, 
is about twice as large as what is observed. But the large shear that is predicted 
would presumably be reduced if, as suggested in the last section, the experiment 
were to be made with a release of latent heat as a function of the computed vertical 
velocity. The consequent ascent of the air at the individual saturated rate and de- 
scent at the unsaturated rate would have the same effect as a large overall reduction 
of the static stability. It would decrease the vertical shear required for given growth 
rates, especially in the region of maximum growth rate of the disturbances. Be- 
cause the equilibrium vertical shear results fromthe balance between the horizontal 
gradient of heating and the horizontal flux of heat, which is mainly the eddy flux 
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Figure 2] 
a) The predicted mean meridional circulation, in cm sv) (Phillips, 1956). 
b) The “observed” mean meridional circulation, incms_ 


_of heat by the disturbances themselves, there would be a reduction of the predicted 
vertical shear. 

The predicted surface zonal winds, in this experiment, although correct in sign, 
are also about twice as large in magnitude as the observed surface winds. Very likely 
the predicted winds would be reduced to about the correct magnitude if, instead of 
the stress coefficient Km = 0.2 x 1074 used by Phillips, Km = 0.42 x 10-4 were used, 
this being our empirically determined value of Km in the absence of large-scale to- 
| pography. 

_ The remaining discrepancy between the predicted and observed mean zonal cir- 
culations is the predicted symmetrical development of the high latitude and low 
latitude easterlies, especially the fact that the low latitude easterlies do not extend 
| to high levels, as they do in the real atmosphere. There are several possible reasons 
for this. One is the difference in the geometries of the two systems, rectangular on the 
‘one hand and hemispherical on the other. A second reason might be the fact that 
in the model the flow is set identically equal to zero at the southern wall. (This boun- 
dary condition could be overcome by integrating over both northern and southern 
‘hemispheres simultaneously, assuming symmetry about the equator.) A third reason 
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a) The predicted mean temperature, in °C (Phillips, 1956). 
b) The observed mean temperature, in °C. 


for the predicted development of the low latitude easterlies may be the constant 
Coriolis parameter that is used in the model. A fourth reason, suggested by Phillips, 
is the neglect of the terms in the vorticity equation which lead to transports of zonal 
momentum by the horizontal and vertical components of the mean meridional 
circulation and the vertical eddy transport by the large-scale (computed) vertical 
velocities of the disturbances. An improvement of the model along these or other 
lines is needed to obtain a more realistic prediction of the low latitude easterlies. 
Figure 21 shows the predicted and “observed”? mean meridional circulations. The 


reason we put quotation marks around the word “observed” is that Figure 21b is 
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not based on actual observations of the meridional wind but is derived as that me- 
ridional circulation which will satisfy the momentum requirements of the real 
atmosphere under the assumption that we can neglect all scales of vertical eddy stress- 
es above the 700-mb level (Mintz and Lang, 1955.)* Having seen that there is fairly 
good agreement (in pattern and in order of magnitude) between the predicted and 
observed zonal surface circulations, and hence between the predicted and observed 
convergences of the poleward momentum flux, we expect to find agreement between 
the predicted and ‘‘observed”’ meridional circulations. Clearly, in the model as in 
the real atmosphere, if we assume no frictional forces in the upper levels, then wher- 


| ever there is divergence in the poleward momentum flux it must be balanced bya 


_ positive Coriolis torque (v> 0), and wherever there is convergence in the poleward 


_ momentum flux it must be balanced by a negative Coriolis torque ‘ < 0). 


In Figures 21a and b, the numerical values show the pressure-averaged values 


_ of v for each poleward and equatorward branch of the meridional circulation. We 


note that in the real atmosphere the intensity of the low latitude direct cell is greater 


_ than the middle latitude indirect one. This is true because the polar cell plays a neg- 


ligible role in the real atmosphere and because the Coriolis parameter is small in the 
tropical latitudes. (As the equations show, the small Coriolis parameter is more 
important than the greater length of the latitude circles. And integrated from 


- equator to pole the total Coriolis torque must be zero: that is to say, the total trans- 
port of momentum from the upper part to the lower part of the atmosphere by 
_ these circulations must be zero.) 


In the model, however, the Coriolis parameter is the same in all latitudes. In addi- 


tion, because of the geometry of the system, the contribution of the polar cell to the 


_ upward transport of momentum (or one can say, its contribution to the total Coriolis 


torque) is not negligible, but equal in magnitude to the tropical cell. Both of these 
reasons operate together to make the intensity of the tropical direct cell smaller than 


that of the indirect cell of the predicted meridional circulation. 


On the right hand side of Figure 21 we show the profiles of ¥ as a function of 


pressure for the central latitude of the tropical direct cells. For the predicted cir- 


culation, this was obtained from the predicted 71 = Vs, given by Phillips, and the 


* We may say in support of this usage that where we do have a reliable determination of the 
mean meridional circulation based on the actual observed winds—in the centre of the tropical 
direct cell and in the winter season—it is very much like the circulation given by the technique used 
to construct Figure 21b. The best of the direct determinations is that of Vuorela (1955), and it shows 


almost exactly the distribution of > as that given for the same latitude and season by Mintz and 
Lang (1955). bi) 
The difficulty with us 1g the direct technique to obtain a reliable estimate of y in the middle and 


higher latitudes is that y is ageostrophic; and in these latitudes, unlike the tropics, the ageostrophic 
component of the observed wind is masked by its much larger geostrophic component. 
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assumption about the distribution of vertical velocity given in our equation (52). 
As we have just indicated, the good agreement between the predicted and observed 
magnitudes of 7 in these tropical direct cells must be attributed to a fortuitous can- 


cellation of errors. 
With respect to the temperature field, Phillips shows only the predicted tempera- 


ture for the 500-mb level. It is a simple matter, however, using the given lapse rate 


(equation 53), to obtain the complete temperature distribution between p; and p3. 
In addition, from (96) and (95), we have 


Ty = S/b + To, 
S=Q+ Qo. 


and 


Using the value of Q(y) given by Phillips, and Qy = 28 x 10° cal m™ day~!, we | 


obtain the values of T4 ~ Ty shown in Figure 22a. 


The predicted temperature field agrees with the observed field (Burdecki, 1955) | 
in its general pattern. But the magnitude of the predicted horizontal temperature | 
gradient, as we have already indicated, is about twice as large as that of the observed | 
field. This discrepancy, as we have said, would probably be eliminated by introducing | 


a term into the heating function such that latent heat would be released when the 
air ascends. This would increase the growth rate of the disturbances, increase 
the poleward heat flux, and reduce the poleward temperature gradient. 


Although there are details that are wrong, the overall remarkable success achieved | 
by Phillips in using the hydrodynamical equations to predict the mean zonal and | 
mean meridional circulations of the atmosphere must be considered one of the land- | 


marks of meteorology. 


2. The mean zonal and mean meridional circulations on an ocean-covered planet, 
given T,(p), at t =0 


In our second (and also third) class of experiments we consider the atmosphere to 


rest upon an earth which is entirely ocean-covered, this ocean having an infinite 


heat capacity and negligible temperature advection so that the ocean surface tempe- 
ratures remain constant with time. 


For computational purposes we simply maintain each Ty at its initial value. But 


we can better understand the physical character of these classes of experiments by | 


substituting in 
Q = b[AT(2,9,1) — (AT)m] 
the expressions 
(AT) nm = Tg,m ier’ T2,m 
and 
AT( 1,91) = T((1,¢) Fe T,(A,9,1) a T(A,¢) ae. [Tom +: T7'(A,¢,1)], (99) 
where T,' is the deviation of T, from its mean global value. 
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Then 
Q = D[T,(A,~) — Ty.m] — bT7'(A,9,t). (100) 


The first term on the right in (100) is a geographically specified forced heating 
and the second term is a temperature dependent heating. 

As a special case we can take the initial ocean-surface temperature as a function 
only of latitude. Then we have 


Q = d(T) — Ty,m] — bT2'(A,9,1). | 


Here the forced heating depends only on latitude and the temperature dependent 
heating varies linearly with the temperature deviations at the level pp. 
A heating function of precisely this kind (although arrived at by following a quite 
| different line of physical reasoning) has been used by Phillips in a second numerical 
| general circulation experiment. The writer was privileged to be shown some of the 
first results of this experiment (made on the IBM 704 Computer at the Massachusetts 
Institute of Technology), but the general and final results have not yet been published 
by Professor Phillips. 


(101) 


3. The general circulation on an ocean-covered planet, given T,(A, y), at t =0 


| In the first two classes of experiment described above there is nothing in the heating 
functions to distinguish one longitude from another. Therefore (unless one introduces 
a large-scale orography) there is nothing to cause time-averaged variations with 
longitude, and any time-averaged velocity field which develops must be zonally 
symmetrical. The first two heating functions, therefore, while adequate for the 
‘study of mean zonal and mean meridional circulations of the atmosphere, cannot 
-be used for the study of the general circulation proper. 

There is no reason, however, why the ocean surface temperature, T,, should not 
be made a function of longitude as well as latitude so that we can use equation (100) 
‘in its general form, 


Q = b[T (A, y) — Tym] — bT2'(A, 9, t)- (100) 


Each column of the atmosphere is aware that it is in rotation about the local 
ertical, but it has no way of knowing which way is north and which way is east. 
‘The west-east gradients of surface temperature and, as far as we can tell, the west- 
‘east gradients of atmospheric heating (cf. Figures 16 and 17) are as large if not 
larger than the south-north gradients. There is every reason to expect that the model 
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atmosphere will respond to a forced heating which is a function of longitude by 
developing a time-averaged velocity field which is not zonally symmetrical but is 
a true general circulation. That is, it will predict what the old terminology called 
“centers of action’. And, of course, the prediction of a general circulation will 
implicitly contain the prediction of mean zonal and mean meridional circulations. 

A numerical general circulation experiment of this kind, using an ocean-covered 
planet with the initial surface temperatures varying with longitude as well as latitude, 
is being carried out by Mr. Avraham Huss, of The Hebrew University of Jerusalem. 
He is using the WEIZAC Computer, at Rehovot, Israel, for his computations. 
The results of this experiment should be most interesting. 


4. The general circulation on a planet covered partly by land and partly by ocean 


Finally, we can carry out a general circulation experiment with a planet partly land- 
and partly ocean-covered. 

Over the oceans, in this experiment, the surface temperatures remain constant 
and equal to their initial values. Here we can take either 77 = 7,(y) or Ty = 7,,(A, ¢) 
so that the heating is given by either (101) or (100). 

Over the land the heating is given by (97), or 


Q = SY, ©) — Qo, 
where © is the time of the year. 

The writer will carry out an experiment of this kind in the near future on the 
IBM 709 Computer at the University of California, Los Angeles. 

Because S varies with ©, the time-averaged motions that will be of greatest interest 
to us in this experiment will be those of the transient general circulation. But 
by averaging over many (annual) cycles, we can obtain a prediction of the normal 
general circulation. 

Assuming that the experiment is successful, it is not at all obvious that in the same 
phase of each annual cycle the predicted transient general circulations will be ap- 
proximately the same. This is one of the things we want to learn from this experiment. 
It may give us some clue to the answer to one of the important questions in meteoro- 
logy, which is whether the large differences that are observed in the seasonal circulation 
(and seasonal weather) from year to year in the real atmosphere are due to extra- 
terrestrial causes or have their origin in the earth-atmosphere system itself. 

It is clear, of course, that our experiment may be deficient in the important respect 
that the ocean, with its great capacity to store heat from one season to another, 
and even from one year to another, is treated in our model as being inert—in the 
sense that in our model it is only the atmosphere itself which determines how much 
heat it will take up from the oceans. We have also neglected other effects which may 
be important in causing the transient circulation to be differeat in different years; 
for example, the variation of albedo, and hence of S, due to variations of cloud- 
or snow-cover over the land. 
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But assuming that after developing a hierarchy of models we eventually find one 
which takes all important effects into account, and gives us in fact a transient 
general circulation that varies from one cycle to another even with a constant external 
energy source, the question will remain whether we can use that model to predict, 
from an observed initial state, the transient general circulation of the real atmosphere. 

The answer to this question will depend on whether, in such a self-contained system, 
the heating which generates the transient general circulation depends upon only the 
broad scale (and therefore observable) features of the initial state, or whether the 
heating depends also upon the details (and therefore for all practical purposes the 
unobservable features) of the initial state. Because the atmosphere is slow in respond- 
ing to its heating, there is reason to be optimistic and to expect that for a forecast 
over some large part of an annual period (say a three-month or six-month forecast 
of the transient general circulation—and also for the forecast of other statistical 
properties of the velocity field, such as the variability of the wind, the cyclone tracks, 
and so on) errors in the details of the observed initial state will be of little importance. 
We can say that with respect to the initial information that is needed (as well as with 
respect to what is predicted), the transient general circulation bears the same rela- 
tion to the details of the instantaneous motion as the temperature of a gas bears to 
the motion of its individual molecules. 

A rational, physical method for long-period forecasts of the transient general 
circulation—that is, for long-range seasonal weather forecasts—is an ambition or 
dream for the future. But we will never arrive at this goal (or even know whether 
we can arrive) unless we begin with the first steps. 
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STUDIES OF THE EFFECT OF SALINE IRRIGATION WATER ON CALCA- 
REOUS SOILS. Il. THE BEHAVIOUR OF CALCIUM CARBONATE 


DAN H. YAALON 
The Hebrew University of Jerusalem 


ABSTRACT 


The presence of CaCO; has often been considered to be a sufficient guarantee for prevention 
of sodium adsorption by clay from saline irrigation water. A scrutiny of the behaviour 
of CaCO; at conditions prevailing in the soil environment was undertaken to test this 
concept. 

It is shown that the physico-chemical equilibria of the soil-CaCO 3 system are analogous 
to the equilibria of the pure CaCO; system. Because of its effect on the pH of the saturated 
solution, there too carbon dioxide pressure of the gaseous phase is the most important 
variable affecting the solubility of CaCO 3. Exchangable ions of the clay exert an influence 
similar to that of dissociated salts. 

It is concluded that under conditions prevailing in the soil during percolation, the solubi- 
lity of CaCO; is too low to supply sufficient Cat +-ions to the soil solution and thus 
unable to prevent the adsorption of sodium. Only in the case of very dilute saline solutions 
is solid CaCO; able to supply sufficient Ca* *-ions to the percolating water in order to alter 
the undesirable Nat+/Cat*+ ratio. The effect of flow rates on dissolution is considered. 

Since gypsum, on the other hand, is sufficiently soluble to prevent the adsorption of 
sodium from moderately saline waters, it is suggested that it be applied as a preventive 
soil conditioner on all soil where permanent irrigation with more or less saline water is 

. established. In order to prevent the dispersion of the clay and to keep the highest possible 
concentration of calcium ions during leaching with salt-free water, it is suggested that 
gypsum be supplied regularly before the winter rains. 


INTRODUCTION 


During percolation of saline water through soils, exchange reactions take place 
between the exchangeable ions of the solid clay phase and the electiolytes of the 
liquid phase, until an equilibrium between the solid and liquid phases is established. 
If slightly soluble salts are present in the soil, they perpetually supply to the percolating 
soil solution additional ions and thus continuously affect the adsorption equilibrium. 
The exchange reactions taking place in calcareous soils during percolation with 
saline irrigation water belong to this special case and will be discussed in the present 
paper. 

The questions asked are: What amount of calcium carbonate can the soil solution 
dissolve, and is this sufficient to prevent the adsorption of sodium by the clay? How 
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does it affect soil pH and other soil properties? What practical conclusions can be | 


drawn? The answer must necessarily be based on the elucidation of the solubility 


relationships of calcium carbonate under the various environmental conditions of | 
calcareous soils, and an integration of its capacity to supply calcium ions to the | 


soil solution with the laws of ion exchange. 


SOLUBILITY AND INTERACTION OF CALCIUM CARBONATE 
IN THE SOIL ENVIRONMENT 


Applying the law of mass action, various workers have shown that in aqueous | 


solutions in contact with solid CaCO; definite relationships are established between 
atmospheric CO, pressure, pH, and ionic strength of the solution (e.g. Bjerrum 
and Gjaldbaek 1919). Based on a systematic theoretical analysis and incorporating 
available experimental data, several equations and equilibrium constants have 


been derived, showing the effect of all relevant variables on the solubility (Yaalon — 


1954 and 1957). 


The system CaCO3-CO2-H,0 at different pH values 


In an invariant system the solubility of CaCO increases logarithmically with decreasing pH and 
is proportional to the cube root of CO, pressure. By combining the effect of these variables and 
the effect of the ionic strength, an overall equilibrium equation (at 25°C) is obtained, 


2pH + log [Catt] + log Poo, = 9.85 + f(u) (1) 


where [Ca* +] is the molar concentration of calcium and f() is some function of the ionic strength 
of the solution. A satisfactory approximation of the theoretical Debye-Hueckel function is given by 
the expression 1.7 +/, derived graphically by plotting ionic strength against the apparent equilibrium 
constant. 

The relation is presented in the solubility curve of Figure 1, demonstrating that CaCO; solubility 
and pH of the resulting solution are determined by the CO, pressure of the gaseous phase. At the 
normal CO2 content of atmospheric air the solubility is only 0.6 mM per litre and the resulting 
PH of the solution is 8.31. 
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Figure 1 
Solubility of calcium carbonate as a function of pH and carbon dioxide pressure (25°C) 


Vol. 7G, 1958 D. H. YAALON 117 


Effect of neutral salts. In dilute solutions up to ionic strength of 0.01, increase in solubility is 
for all practical purposes a linear function of /# tesulting from increments of neutral salts. At 
higher concentrations the relation is curvilinear. Generally the increase in solubility is accompanied 
by a slight decrease in the pH of the saturated solution. At the concentration of electrolytes and 
at the carbon dioxide pressure likely to occur in soils, the increase in solubility is slight and will 
not be over 1 mM per litre. 

Effect of a common ion. A salt with a common ion decreases the solubility of CaCO3;, but on the 
other hand, if the common ion is Ca++, it permits a much higher total calcium concentration in the 
solution. This in turn is limited only by the solubility product of the more soluble salt. Hence, if 
CaSO, is also present in the soil the calcium-ion concentration can reach a value of about 30 mEq 
per litre, and is practically unlimited if the salt is CaCl». 

When considering the interaction of calcium carbonate with the clay minerals of a soil, the impor- 
tant variables affecting it are (1) the carbon dioxide content of the soil atmosphere, (2) the ionic 
composition and strength of the soil solution, and (3) the exchangeable ions of the clay complex. 
The effect of these variables will be examined in the following. 


The effect of carbon dioxide 


The effects of variations in carbon dioxide pressure on CaCO; solubility are large in comparison 
with other factors.. 

Though at constant pH an increase in carbon dioxide pressure causes a decrease in solubility 
(cf. equation 1), when pH is permitted to adjust itself, higher pressure of CO, is accompanied by 
a decrease in pH (cf. equation 3 below), with an overall effect of enhanced CaCO, solubility. Thus 


log [Cat +] = 1/3 log Pco, — 1.95 +f(u) (2) 


The constants of equations 1 and 2 are the product of the equilibrium constants of the participating 
simultaneous equilibria, and have been calculated for 25°C (Yaalon 1954), using the best available 
thermodynamic constants. Figure 1 shows graphically also this relationship. 

Consequently, the carbon dioxide content of the soil atmosphere must be regarded as one of the 
important variables affecting the capacity of CaCO; to release Cat *-ions to the soil solution. The 
range of its concentration in the soil atmosphere is a function of biological activity and is therefore 
much wider and higher than that of the air above the soil. It is known to vary between the concent- 
ration of the atmospheric air, 0.03%, and up to 5-10% (Russel 1950, p. 339). Not many exact 
determinations of CO, content in irrigated arid soils are available, but assuming it to be somewhat 
lower than in humid soils, its normal range can be tentatively set to be between 0.2-1.0%, higher 
values being applicable for tropical soils with considerable microbiological activity. At this range 
of CO, pressure, the solubility of CaCO; reaches 1.0-1.9 mM per litre and imparts to the solution 
a pH of 7.3-7.8, which is the normal range of pH in calcareous soils under field conditions (Figure 1). 

The pH of a CaCO; solution in equilibrium with atmospheric air (Pco, = 0.00032 atm) is 8.31 
at 25°C, which should therefore also be the pH of calcareous soils aerated with atmospheric air 
of normal composition. Usually it is found to be somewhat higher, but does not exceed 8.5 (Bradfield 
and Allison 1933). This might be due to the inherent colloidal properties of the clay, but is generally 
ascribed to the presence of exchangeable monovalent bases, usually sodium. Also MgCoO;3, which 
is more soluble than CaCO; and yields a higher pH at similar CO, pressure (Gjaldbaek 1921), will 
contribute to this effect. 

The effect of carbon dioxide on the soil reaction of calcareous soils has been investigated by 
Whitney and Gardner (1943). The curve showing the effect of partial CO, pressure on pH was found 
to follow the same logarithmic shape as for pure carbonates. The average regression for the 18 soils 
tested has been calculated by the present author and results in the following equation: 


pH = —0.65 log Pco, + 6.00°(25°C) (3) 


i i i i i Icium carbonate, which are 
The constants are practically identical with those obtained for pure ca ate, 
9/3 and 5.90+f WD respectively. A more careful scrutiny of Whitney and Gardner’s data reveals 

i iff i i i lon 1957). 

lightly different functions for three different types of soils (Yaa 
: 78 in addition to evidence obtained recently by Turner and Clark (1956) and Cole (1957), 
Sriables us to conclude that the effect of carbon dioxide on CaCO3 solubility in calcareous soils is 
analogous to its effect on the pure salt, and can be deduced with sufficient accuracy from theory. 
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The effect of the soil solution 


The concentration of the soil solution, defined as the solution existing at a soil moisture content 
within the field range and containing the salts of unadsorbed anions, varies inversely with the amount 
of moisture present in the soil. The change in total concentration is accompanied by a shift in the 
relative concentration and distribution of the different ions, mainly as a result of different absorbing 
powers of the exchangeable cations in the ionic double layer of the soil phase (Reitemeier 1946). 
Consequently there exists no static boundary between the adsorbed ions in the solid phase and those 
in the soil solution. 

The soil solution of a leached agricultural soil under normal conditions of moisture contains up 
to 500 ppm of soluble matter (Russe! 1950, p. 438). Reckoning with the normal composition of the 
common ions, this becomes equivalent to about 7 mEq/litre with an ionic strength of « = 0.01 and 
osmotic pressure of 0.25 atm. Hence the soil solution in such a soil is very dilute, even at low moisture 
content. 

In soils of arid regions, because of the occurrence of slightly soluble salts and minerals and because 
of retarded leaching, the concentration is usually several times higher. In calcareous soils the concent- 
ration is about 1000-5000 ppm and can be much higher in calcareous-saline soils. Reckoning with 
the above range, we obtain an ionic strength of ~ = 0.04-0.15 (4/ = 0.2-0.4) (expressed in equi- 
valents 20-100 mEq/litre), and an osmotic concentration of 1-4 atm. The solubility of CaCO, at this 
range of ionic strength and normal CO) pressure is 4-5 mEq/l. and the pH of its solution 7.2-7.7. 
But if due to the presence of Nat-clay the pH of the soil rises above 8.0, and carbon dioxide pressure 
is decreased as a result of reduced biological activity and air exchange, then the solubility of CaCO; is 
less than 2 mEq/1. 

The analytical results of Ca* +t concentration in the soil solution of the permeability experiments 
reported in Part I of this study (Yaalon 1955) agree well with the above considerations, as do also 
concentrations of soil solutions at field moisture and saturation percentage as determined by Reite- 
meier (1946). For calcareous soils not containing gypsum, the concentration of soluble calcium 
was within the limits discussed above. Where Na2CO; was present, and the pH was high, the concent- 
ration of calcium was considerably lower, while CaSO, was able to keep the concentration of calcium 
at alevel of about 25 mEq/litre. 


The effect of exchangeable ions 


The most important variable affecting the behaviour of CaCO; in the soil is the exchange complex 
of its clay minerals. The calcium dissolved from the calcium carbonate may be exchanged for other 
cations: 

Cat+ + M-clay =— Ca-clay + M++ (4) 


Thus the solubility equilibrium of the CaCO; -soil solution system (equation 1 and 2) is affected 
by the additional equilibrium with the clay mineral (equation 4). Hydrogen, calcium and sodium can 
typify the main types of exchangeable ions encountered in the soil, and their effect on the behaviour 
of CaCO; will be considered. 

If hydrogen is the exchangeable ion of the clay, the acid soil will cause the decomposition of the 
CaCO; with a resulting neutralisation of the clay, a reaction utilised in the liming of acid soils in 
humid areas. The solid of CaCO; constitutes a reserve of calcium ions, which will tend to keep the clay 
saturated with bases and maintain a constant pH. The level of the pH will depend solely on the 
carbon dioxide pressure and on the electrolyte concentration of the soil solution. Only when all 
CaCO; has been dissolved will there exist the possibility that the pH decreases, since the buffer 
capacity of the bicarbonate system is rather small. ; 

Cat *-clay willinteract in a manner similar to that of an ordinary ionised Cat++-salt at the same 
Cat t-ion activity. A “titration” curve of CaCO; with increasing amounts of Cat +-clay is presented 
in Figure 2. The data were obtained at constant Pco, of 0.00032 atm (by aeration) and at 25°C. It 
can be seen that the effect of exchangeable Catt on the pH and solubility is very small which is 
in agreement with Marshall’s determinations of the active fraction (Marshall 1950). For ts Cat +. 
clay an activity coefficient of not more than 0.09 was obtained, kaolinite having the highest value 

As exchange reactions between Na+ and Cat * will take place when Na*-clay is added to a CaCO 
solution, it cannot be a priori expected that the Nat-clay will behave similarly to a disocinthd 
neutral Na*-salt. Because of its alkaline reaction, sodium clay will also effect solubility differently 
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than a neutral Na*-salt. A “titration” of CaCO; with Nat-clay, analogous to Figure 2, is presented 
in Figure ae ot can be seen that increasing amounts of Nat*-clay increase the PH of the resulting 
equilibrium solution, the extent of the increase depending on the CO, pressure. Also, in agreement 


with the effect of pH on solubility of the pure salt, the calcium ion concentration in the solution 
decreased. Similar curves have been obtained by Puri and Puri (1938). 
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pH of a saturated CaCO; solution as a 
function of added Ca-clay (at CO, 
pressure of atmospheric air) 


PH of asaturated CaCO; solution as a 
function of added Na-clay (at CO, 
pressure of atmospheric air) 


In conclusion, it can be said that the physico-chemical equilibria of the CaCO; system are not 
appreciably affected by the presence of soil minerals, and hence that the laws governing these equilib- 
ria are directly applicable to the soil CaCO system. Collaterally it is possible to interpret observed 
soil phenomena in terms of solubility, pH and ionic strength. 


THE BEHAVIOUR OF CALCIUM CARBONATE DURING PERCOLATION 

The preceding chapter has demonstrated that under the environmental conditions 
prevalent in calcareous soils the solubility of CaCO; is low, and that even under 
relatively favourable conditions it would reach only about 5 mEq Ca*™ per litre. 
CaCO; is thus ordinarily unable to supply sufficient calcium ions to a saline solution 
to alter materially an unfavourable proportion of Na‘/Ca**. Normal exchange 
reactions between the percolating liquid and the exchangeable ion will lead to a 
giadual saturation of the exchange complex with sodium, which may eventually 
attain an undesirable level of saturation (Yaalon 1955). 

Having derived the amount of CaCO; carried in the soil solution in equilibrium 
with a given environment, there remains to consider the rate at which this equilib- 
rium concentration is attained. 

The rate of dissolution of a solid is related to its solubility ; the greater the solubility, 
the higher the rate of solution. Of the many variables affecting and controlling the 
rate, the more important ones are: (1) the exposed surface and structure of the solid, 
(2) the rate of diffusion of the dissolving ions, which in turn is a function of the concent- 
ration gradient and theionic mobility, and (3) the state of motion of the dissolving liquid. 

Consider a granule of CaCO; embedded in a column of soil through which water 
is percolated. The percolating liquid removes dissolved ions from the vicinity of 
the granule as long as a saturation deficit exists. The maximum rate of dissolution 
is encountered at. high rates of flow with solutions possessing a high capacity to 
hold calcium. When the percolation rate decreases—in soils primarily due to the 
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swelling and dispersion of clay—the rate of dissolution decreases for two reasons. 
Asthe rate of removal of the solution film surrounding the calcium carbonate decreases, 
the concentiation gradient is gradually reduced, with a consequent decrease in the 
diffusion rate of the ions. Also, as a result of ion exchange and hydrolysis, the capacity 
to hold calcium is further decreased. 

Permeability characteristics of the soils reported in Part I of this study (Yaalon 
1955) agree well with these considerations. In all cases where permeability has been 
reduced by the swelling of the clay, the resulting slower rate of flow, the water being 
forced to flow through narrower and more twisted channels, has reduced the rate 
of dissolution of CaCO; to such an extent that it has become ineffective as supplier 
of calcium ions to the solution. 

Information on the amount of CaCO; dissolving during percolation with saline 
solutions was obtained from the following experiment. NaCl solutions were percolated 
through a 10 cm deep bed of precipitated CaCO3 (< 1m) and washed quartz sand 
(0.05-—0.imm), thorcughly mixed at the weight ratio of 1:2. The calcium ion 
concentration and pH of the percolate were determined. In Table 1 the concentrations 
reached during steady flow are compared to the concentration of the perco ating 
solution used. Considering these values as the limiting values of Cat * concentration 
in the soil solution, it is evident that only where the percolating solution is very weak 
will sufficient calcium carbonate dissolve to alter its undesirable ratio of N*/Ca‘*. 
A solution containing 6 mEq Na* and 3 mEq Ca‘* per litre (67% Na*) was 
percolated through the CaCO3 sand bed and passed through practically unaitered. 


TABLE I 
Effect of percoiating solution on composition of leachate 


Medium: Precipitated CaCO; mixed with quartz sand (0.05-0.1 mm) at weight ratio 1:2 
Percolating fluid: NaCl soiution. Room temperature and Pco, of atmospheric air 


Concentration of Composition of leachate at steady flow conditions Equivalent Nat 
percolating solution percentage 
NaCl Cat+ pH of leachate* 
(mEq/1.) (mEq/I.) 
0 0.85 8.70 _ 
1 1.18 8.45 46 
D 1.03 8.65 66 
5 0.97 8.90 84 
10 0.98 9.10 91 


* Tonic strength (./u) of leachate is about 0.12-0.16 


PRACTICAL IMPLICATIONS 
To prevent the adsorption of sodium from irrigation water with an undesirable 
Nat/Catt ratio, a considerable increase ir the calcium ion concertiation of the soil 
solution must be obtained. The logical conclusion for calcareous soils would be to 
try to increase the solubility of the lime present, which after all presents an almost 
inexhaustible resource of calcium. 
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Since carbon dioxide is the most important factor enhancing lime solubility, soil 
management and irrigation practice must aim to fulfill this requisite. Improved 
aeration in itself will not increase the carbon dioxide content of the soil atmosphere 
unless it is accompanied by increased microbiological activity or tespiration of the 
root system. Hence plants with an extensive root system are preferred, and a periodic 
addition of organic manure should be beneficial. Alsc, the evidence that drying 
of the soil improves permeability should be examined in this light. 


Should the above measures not suffice, and a decrease in permeability be observed, then means 
of increasing the calcium-ion concentration by a suitable external addition must be considered. 

Gypsum is many times more soluble than CaCO; and is able tosupply about 25-30 mEq Cat + ions 
per litre. Its capacity to alter the composition of the liquid phase is therefore considerable. Its solubi- 
lity relationships in soils have been studied by Vanoni and Conrad (1942) and found to be explicable 
in terms of solubility product and ionic strength phenomena. It fulfills all requirements of a suitable 
amendment, and can be added either to the irrigation water or directly to the soil. 

It seems to the writer that in arid and semi-arid irrigated agriculture, gypsum application should 
take the place of liming in humid climates, and that it should be introduced as a regular feature 
in soil management, especially where a permanent system of irrigation agriculture has been established 
and water of more or less doubtful quality is being used. 

It is probable that better results will be obtained with smaller and more frequent additions. An 
application prior to a change in the composition of the water—this includes the effect of salt-free 
winter rains—should be beneficial and enable the establishment of a calcium ion concentration at a 
sufficiently high level to keep the soil colloids in a flocculated and permeable condition. 
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SIMULATION, BY DYNAMICAL MODEL, OF SAND TRACT 
MORPHOLOGIES OCCURRING IN ISRAEL 


M. Rim 
The Hebrew University of Jerusalem 


ABSTRACT 


The question of which configuration is attained (and maintained) by the surface of a bed of 
loose sediment when it is subjected to the erosive and simultaneously depositional action 
of a sediment-carrying stream is studied with the aid of a hydraulic flume. The experimental 
results obtained for this particular system have significance also for analogous systems 
(where the fluid may be a fluid other than water) if a certain relationship exists between the 
scales for geometrical dimensions, velocities, forces, specific gravities, viscosities, etc., in the 
experimental system (model specimen) and the corresponding scales in the analogous 
system (prototype specimen). This relationship between the scales is defined here in a 
non-rigorous manner (by three “Prerequisites of Dynamical Analogy”) which makes 
it applicable to geomorphological and geological problems. 

On the basis of these Prerequisites the flume experiments and the miniature landscapes 
created in the flume under the water current are postulated to be dynamically analogous 
to large-scale landscapes formed by wind—which analogy is then proven to be real by 
actual comparisons. Thus it becomes possible to study the mechanism of geomorphological 
evolution and the relative importance and interdependence of the factors involved on a 
dynamical model instead of on a full-size dynamical prototype, which latter would usually 
be impracticable. Here such a study of a collection of dynamical prototypes, the typical 
sand tract morphologies occurring in Israel, is reported. 

Dynamical model studies of fossil landscapes could also contribute to the identification 
of wind regimes and environmental conditions of past periods. 


THE PREREQUISITES OF DYNAMICAL ANALOGY 


Dynamical models are functioning models in which, besides geometrical proportions, 
the time scales and mass scales (and thereby the inertia and force units) also bear 
a fixed relationship to the imitated system, which is then referred to-as the prototype. 
Designers of machinery and structures involving relative motions of fluids and solids 
(such as missiles, airplanes, ships, river channels, wave- and wind-breakers, etc.) 
depend fully on dynamical modelling, and it is also a dependable tool in geophysical 
- investigations that deal with problems concerning atmospheric or oceanic currents. 
As long as the fluid involved is of the “Newtonian fluid” kind, and as long as the 
confinements are not deformable by the fluid’s motion, the correspondence between 
a dynamical model and its prototype is mathematically rigorously definable by the 


“criteria of dynamical similitude’’.* 


* A set of formulae that are derived by ‘dimensional analysis” of the relative parameters 
and which shall make it possible to predict, for instance, whether a half-size model airplane of 
say one-tenth the weight of the full-size prototype will lift itself off the ground at one-half of, 
the same as, or twice the prototype’s runway velocity. 

Received July 18, 1955. 
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Geomorphological prototypes, on the other hand, present mathematical difficulties 
because they do not fulfill the two above-mentioned conditions. Firstly, the fluids 
involved, which in geomorphological prototypes naturally are always carrying 
sediments, no longer display “Newtonian” behaviour. Secondly, geomorphological 
systems deform their confining boundaries by erosion and deposition. 

Thus it is clear that some of the mathematical rigour of the criteria of dynamical 
similitude must be sacrificed when interpreting analogies in geomorphology.* 

For the task confronting the author in the work reported here it was sufficient to 
simplify the rules of dynamical analogy and reformulate them as ‘General Prerequi- 
sites of Dynamical Analogy in Eroding and Depositing Systems” which state the 
conditions that must be fulfilled if the functioning of a prototype is to be re-enacted 
on different time-, size-, and mass-scales. 


Prerequisite 1 


When a bed of loose sediment is exposed to the action of a fluid flowing over it, the 
equilibrium between erosion and redeposition will lead to a steady (not static) state 
of configuration of the bed surface. The height of the relief of the bed surface, usually 
a pattern of undulations and knolls, is then commensurate with the height of what 
will be termed here “‘sediment atmosphere,” namely the churned-up sediment cloud 
that consists of the fluid-carried portion of the sediment (suspended, saltating or 
otherwise moving) and blankets the solid sediment bed. There is a fixed relationship 
between relief and atmosphere heights in dynamically comparable systems; in 
other words, if R,; and R» denote the relief heights in two systems and A; and A) 
correspondingly the atmosphere heights in these systems, then Rj/R2 = A,/Az, 
irrespective of fluid velocities. 


Prerequisite 2 

Since the sediment atmosphere is isotropically turbulent except nearest to the sediment 
bed, it determines the horizontal dimensions, i.e. the dimensions of the plan pro- 
jection of the sediment-bed surficial features, in such a way that A,/A, » L,/L,, 
where L stands for corresponding distances between markers on the bed surfaces. 
For the lateral distances between any other markers, L’, L’,L”,...L™,..., the 
same relationship holds, of course. Consequently, 

ES OS el ARE A Be Ss ag BUR FRae ey oy 
i.e. there ought to be full geometric similarity between dynamically similar systems. 
Conversely, full geometric similarity between the surficial patterns formed in two 
systems is a weighty argument, though not a fully conclusive one, in favour of the 
assumption that these two systems are functioning dynamically analogously. 


“89 


* It is possible to conceive of the confining sediment bed not as deformable walls and floor, 
but as the multitude of rigid grains it really is, and to treat the fluid not asa sediment-carrying non- 
Newtonian mixture but as a clear Newtonian fluid that surrounds each suspended grain. This approach 


leaves the problem in principle a rigorously solvable one, but the solution would be impractically 
complicated. 
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Prerequisite 3 


In the case of non-steady multi-directional fluid flow, the “resultant” flow vector 
does not necessarily coincide with the “resultant”? vector of sediment drift. The 
latter depends on a host of environmental factors, particularly on local differences 
of sediment configuration and abundancy. For two systems to function as dynamical 
analoga it is necessary that the lateral distribution of sediment drift be similar in 
both systems. Conversely, similar drift distribution in two systems may indicate 
analogy of dynamical functioning in other respects. 


Additional Prerequisites and Simplifications 

If the systems contained sediments other than completely loose sand, it would 
have been necessary to add here prerequisites concerning the relative scales of the 
cohesion forces in the compared systems. This problem is too complicated and 
therefore all modelling studies so far have been confined to modelling eolian sand 
tract morphologies. A prior dynamical model study of eolian sand landscapes by 
hydraulic flume (a study of the Sahara desert sand dune morphologies) had been 
completed with the aid of similar experimental aids by L. Escande (1950, 1951). 
Escande’s prototypes were selected as the most simple (or assumedly simple) existing 
systems: he took for granted in his prototypes that winds were unidirectional 
and that sand was omnipresent. These assumptions eliminated the problems connected 
with Prerequisite 3 and made it possible for him to apply formally, without modifica- 
tion, the mathematical criteria based on dimensional analysis. By means of these 
criteria Escande was the first to show that sub-aerial sand landscapes must be model- 
led dynamically not by miniature sand landscapes formed in wind tunnels, but by 
those formed under flowing water. This is what is expressed implicitly in the first 
of the above prerequisites. 


METHOD AND RESULTS OF EXPERIMENTATION 


The hydraulic flume in which the dynamical model experiments were carried out 
consisted of a 10 m long channel of square cross section (0.35 x 0.35 m) and a pipe 
fitted underneath through which the fluid was continuously pumped from the channel 
tail end back to the channel head (Figure 1).* The channel fioor was completely 
covered with fine sand, forming the sediment bed. The surface of the bed attained and 
maintained specific shapes when water was made to flow through the flume channel. 

The variables were the kind of sediment, sediment distribution and current velocity. 
Only one kind of sediment was used in the series of experiments reported here, 
namely, fine sand consisting of well-rounded quartz grains of uniform size (sieved 
and winnowed by means of a specially constructed machine which removed all 
grains that were larger than 0.14 mm or smaller than 0.11 mm in diameter, or deviated 


* This flume was built for sediment traction studies by the U.S. Soil Conservation Service, at 
the Sedimentation Laboratory of the California Institute of Technology. The author thanks Profes- 
sor V. A. Vanoni for the permission to use this flume, and for the help and instruction received 
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too much from spherical shape or average specific gravity). Local sediment distribution 
in the flume was influenced (in a manner described in detail further on) by placing 
“false bottoms” into the channel. Water-current velocity was regulated to the desired 
values by inclining the whole channel lengthwise to the horizontal (and adjusting 
correspondingly also the pump speed in the re-circulation pipe). The maximum 
useful velocity in this system was 1.2 m sec as measured in the channel centre, beyond 
which critical velocity the current attained the undesirable characteristics of “shoot- 
ing flow’*. In the following, the expression ‘“‘swift current’’ designates all currents 
of velocities between critical velocity and half this velocity. “Medium current” 
designates all currents of velocities of between one-half and one-fourth of the critical 
velocity. “Slow current” refers to all currents of velocities below one-fourth of the 
critical one and above threshold velocity (the velocity at which the current starts 
to be able to move grains at all). The latter velocity depended very much on the 
momentary bed smoothness, and was about one-tenth of the critical velocity. 


Experiment Series 1: Effect of swift current 


A swift current allowed to act on the sediment bed lifts up so much sediment that 
the total channel height of fluid is turbid with sediment and always A/R > 1. 
It establishes in the course of only a few minutes a bed configuration pattern that 
is in harmony with the eroding and depositing action of the current: straight, trans- 
versally aligned equidistant corrugation (see Figure 2, where these are seen as photog- 
raphed through the glass side-windows of the flume). Preferential erosion takes place 
on those slopes of the corrugations that face away from the current and preferential 
deposition on those that face into the current, so that the corrugations, while main- 
taining their shape, translate themselves slowly in the same direction as the current 
flows, at about one-tenth the current velocity. The height of the relief increases some- 
what with current velocity. 

No experiment, either in this series or in any of the others, ever made it appear 
likely that there exists a set of conditions where the equilibrium relief could be made 
nil, i.e. where a flat bed would remain a flat bed. (Even shooting flow does not sweep 
the bed surface into a flat shape, but digs into it at some unpredictable spots and 


builds it up at others.) Theory too negates the possibility of dynamical SOY of 
a flat bed. 


Experiment Series 2: Effect of medium current 


Medium currents raise loose sediment from the bed just enough to make the lower 
portion of the flowing fluid turbid, so that A/R w 1. 


* So-called shooting flow ensues when current velocity reaches the order of magnitude of gravity- 
wave propagation. Then the moving fluid elements have no chance of adjusting their motions into 
a geometrically repetitive pattern. (Normally such a pattern is effected by gravity waves running ahead, 
being reflected from obstacles and running back. Under shooting-flow conditions there is only 
one-way Communication between the bed configurations on different parts of the channel and no 
true equilibrium pattern results.) 
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To attain a bed configuration that persists (i.e. perpetuates itself instead of giving 
way to new evolving patterns one after the other), up to several hours of recirculating 
the current over the bed may be necessary, if the current chosen is a medium current. 
The self-perpetuating pattern which eventually emerges does not show that high 
degree of simplicity and uniformity displayed by the patterns formed by swift current. 
It is conceivable that there exist some less simple starting conditions of bed shape 
that lead into a different cycle of evolution than the one obtained in the author’s 
experiments described in this paragraph (e.g. see Experiment Series 6 further on). 
The simplest starting conditions are a bed shape in equilibrium with a swift current, 
the speed of which is then reduced to the lesser speed of a medium current. Then the 
troughs of the transversal pre-existing corrugations are seen to become gradually 
filled up with sediment, but only at equi-distant nodule points, about one corrugation 
distance.apart. The end product observed was an array of rectangularly arranged 
erosion holes as shown in Figure 3. 

This array translates itself in the direction of the current at a small fraction of the 
current speed by the following observed mechanism. Suddenly, at irregular time 
intervals and at unpredictable spots, protuberances of suspended sediment are seen to 
shoot out of an erosion hole and then dissipate, while being carried down the current, 
by diffusing and settling into the sediment atmosphere. These protuberances are the 
result of spontaneous local increase of current turbulence. On closer inspection it is 
seen that the part of the erosion holes most afflicted are the slopes facing the current. 
These slopes are thus-eaten up by erosion exceeding deposition while the opposite 
slopes are built forward by the complementary depositional activity of the current. 

A variant of the above experiment series consisted of using non-steady medium 
current, such as intermittent swift current. The. resulting bed shapes turned out 
to be similar to the above-described ones, often even of improved uniformity and 
regularity. 

Experiment Series 3: Effect. of slow current 

Slow currents do not create a sizeable sediment atmosphere, but merely propel the 
surficial sediment grains rolling or sliding over the bed. Using the terminology of the 
sediment atmosphere, this means that A/R < 1. Under these conditions the initial 
abundance and local distribution of sediment on the bed surface becomes a decisive 
factor, so that one cannot give a generally valid description of the effect of slow 
currents. Even for the case of omnipresence of sediment on the channel floor and 
complete symmetry of the current, one cannot say for sure that evolution of the bed 
surface shape will lead to a self-perpetuating pattern of defined size. 

When slowing a medium current that flows over a bed it has shaped (into the 
shape described in the preceding paragraph) the following changes can be observed. 
The crests bend gradually around so that the erosion holes become circular, and 
at the same time dwindle away at the down-current side of each erosion hole, so that 
the crest system finally breaks up into individual semi-circular knolls, each knoll in 
front of a flattened-out erosion hole. These knolls grow slowly in a haphazard manner, 
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the larger ones at the expense of their smaller neighbours. Their shape, as seen from 
the top, retains more or less the shape of a crescent moon or a horseshoe. Translation 
of the pattern units in the current direction is individual; the more the features grow, 
the slower the translation. There seems to be no definite end to this growth but 
gradual freezing of all translational motion of the surface features with increasing 
relief. Small-size features (ripples) that originate spontaneously here and there 
translate themselves with considerable speed over the current-facing slopes of the 
large knolls (see Figure 4). 


Experiment Series 4: Effect of sediment scarcity 

If there is a limit to the availability of sediment on the channel bottom so that a non- 
erosible floor under the sediment is exposed in the low spots of the relief, knolls 
cannot grow without limit at the expense of their neighbours since they become 
isolated by the baring of the low portions of the relief. Evolution continues then by 
virtue of the following mechanism. As small knolls travel forward faster than big 
ones, the small ones continuously overtake the bigger ones. Whenever it happens. 
that a big one is in the way of a small one the two coalesce. The final self-perpetuating 
pattern is reached when all knolls have identical size and thus identical travel speeds, 
which state is reached in a surprisingly quick time. In Figure 5a, an almost final 
stage of an array produced* in this way of horseshoe knolls on a non-erosible 
substratum is seen photographed from the top through the water current. 


Experiment Series 5: Effect on non-symmetry of sediment availability 


If during an experimental run, sediment was plentiful all the time on one, say the left, 
bank of the main current bed and continuously trapped on the other bank, sediment 
would drift also, energised by the fluid motion’s random components (turbulence- 
agitated diffusion), across the current direction. Such conditions of non-symmetry of 
sediment availability could be created simply by tilting slightly sideways the false 
bottom described in the preceding paragraph. The effect of cross-current drift of 
sediment on the sediment-morphology on the bottom is very pronounced. With 
plentiful sediment and stronger currents, the sediment shapes into obliquely aligned 
corrugations similar to the artificial corrugations shown in Figure 5b. With scarce 
sediment and slower currents, the pattern of horseshoe knolls shown in Figure 5a 
is obtained, but distorted so that the right prongs of the horseshoes are extended 
and the left ones reduced to rudiments. The degree of this distortion is a measure of 


* For convenience, sediment-scarcity in these experiments was not produced by emptying the 
flume and connected pipes of a portion of the sand, but by the following artifice. A “false bottom” 
with legs was placed into the channel so that the sediment-carrying current flowed partly underneath 
it and partly over it. The level of this false bottom was variable and so was the configuration of the 
bottom itself. In the experiment described in this paragraph, the false bottom was flat only slightly 
roughened. In experiments described in the following paragraphs, the false bottom was corrugated 
as shown in Figure 5b. The length of the false bottom was 150 cm (the length of the flume’s working 
section) and its width 8 cm narrower than the flume’s width. Its leg height was adjusted to give the 
desired degree of sediment scarcity in each case. 
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the intensity of the cross-current sediment drift. In Figure 5a some slight distortion 
of this kind can actually be perceived, owing to some extremely slight tilt of the 
false bottom during production of this particular specimen. 


Experiment Series 6: Effect of pre-existent bed configuration 
A false bottom of non-symmetrical non-erosible shape (a rectangular piece of 
corrugated sheet, cut so that the corrugations run 45° inclined to the edges and 
current direction, as shown in Figure Sb) inserted into the current in the flume proved 
to induce cross-current drift of sediment even if sediment availability itself was pains- 
takingly made equal at both sides of the current. The drifting could be noticed by eye. 
In the course of time, however, current transversal ripples superimposed themselves 
on the corrugation crests and on the corrugation troughs. Growing steadily, they 
tended to change the oblique alignment of the corrugations to transversal alignment. 


SAND TRACT MORPHOLOGIES OCCURRING IN ISRAEL 


The sand tracts occurring in Israel, except those scattered among the mountainous 
stone desert in the south of the country, are coastal dune areas of different morpholo- 
gical age, substratum, surroundings and wind regime, displaying different morpho- 
logies in different regions. Their common denominator is that the sand source is a 
defined line, namely the beach. The sand morphology specimens described in the 
following shall be referred to as prototype specimens in view of the fact that the 
experiment series reported here had been selected to model dynamically just those 
typical specimens. The specimens shall be presented here in such sequence that 
Experiment Series 1 models Prototype Specimen 1, Experiment Series 2 models 
Prototype Specimen 2, and so on. 


Prototype Specimen 1: Ripple marks 


The sediment atmosphere for the case of wind-driven sand is determined by saltation 
height and is thus fairly constant, about 30 cm, unless there is very strong wind 
(or local turbulence for other reasons). Thus, eolian morphological features that 
are dynamically analogous to those water-produced morphological features for 
which A/R > 1 (Experiment Series 1) are those that have the smallest geometrical 
dimensions: the eolian ripple marks found on every wind-swept sand surface as 
seen, for instance, in Figures 6a and 6b (which are not typical of any specific locality. 
The white arrows indicating wind direction show the transversal alignment of the 
ripple marks, which typically is also the alignment of their dynamical counterpart 
shown in Figure 2). They are mentioned here only for the sake of completeness. 


Prototype Specimen 2: Erosion hollows 
A region where sand is transported not by steady wind but rather by intermittent gusts 
from opposite directions (the turbulence of which does succeed in raising temporarily 
the sediment atmosphere considerably) is the plain of Acre. A specimen of the sand 
morphology found on this plain is seen on the aerial photograph in Figure 7. A dia- 
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gram showing the wind regime in winter (full lines) and in summer (broken lines) as 
stated in Ashbel’s Meteorological Atlas is reproduced in the marginal inset of the same 
photograph. None of the wind directions is predominant, since the shifting winds 
mutually cancel their direct shifting effects. Nevertheless, sand travels inland steadily, 
owing to the sand-trapping action of the swampy and vegetated land on one side 
and the continuing sand supply from the beach on the other side, a sort of diffusion 
being effected by the random components of the winds. Conditions here resemble 
best the ones produced in the variants of Experiment Series 2, and so do the morpho- 
logical features obtained. There are zig-zag crests enclosing erosion hollows, so 
that the latter form a rectangular pattern of their own, diagonally aligned to the 
general direction of the crests. In winter the wind gusts are seen to whip up sand 
protuberances from the hollows, which presumably are the dynamical prototypes 
(several hundred times larger) of the protuberances seen in the model experiments. 
Also all geometrical dimensions of the prototype features are several hundred times 
larger than those of the model, so that A/R * 1. 


Prototype Specimen 3: Horseshoe dunes 


The wind regimes in the rest of the prototypes to be described here are characterised 
by the dominance of one steady wind direction (see diagram on the marginal inset 
of Figure 8). Thus the sediment atmosphere is always low and sand would travel 
mainly in the dominant wind direction, analogous to conditions of Experiment 
Series 3, if sand were everywhere deep and plentiful. However, there is no large 
continuous area where sand depth would not be a limiting factor, mainly because 
the ground-water table is everywhere in this coastal region close to the dune surface 
and the water-pervaded substratum is not erosible by wind. Thus there is here 
no real continuous large area covered with ever-growing dunes comparable to those 
in Experiment Series 3 produced on model-scale. (Such dunes would look like those 
shown on Figures 6a-and 6b, centre and prong respectively.) Only isolated dunes of 
this kind occur in odd patches here and there. 


Prototype Specimen 4: Array of horse-shoe dunelets 


The non-erosible substratum mentioned above as the growth-limiting factor of 
horseshoe dunes causes the incipient horseshoe dunes to form a repetitive array 
which is persistent, since the individual horseshoes do not have a chance to grow 
at each other’s expense. They may grow only until the substratum becomes exposed 
in the deep spots of the array. To show the array on a single plate, vertical aerial 
photographs have been fitted together in Figure 8 to form a “mosaic”. (The original 
photographs can still be distinguished by the differences of tone on the two sides 
of the straight lines where they were joined). Scale and wind regime are indicated 
diagrammatically in one of the marginal insets. In the other inset the position of the 
sand-source (the coast) is indicated, which extends parallel to the bottom edge of 
the mosaic, about 500 m beyond it. 
This is rather similar to the effect produced in Experiment Series 4. 
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Prototype Specimen 5: Elongated dunes* 

The horseshoes seen in Figure 8 display asymmetry similar to the asymmetry of 
those experimentally obtained and shownin Figure 5a. Their right prongs are more 
developed than their left ones. As in the case of the model, this phenomenon is asso- 
ciated with asymmetry of sand availability. The asymmetry in the prototype case 
is attributable to the sand-supplying coast being on one side of the main wind stream, 
the sand-trapping inland on the other. Thus there is traffic of diffusive sand motion 
across the main wind stream. 

In Figure 9, which depicts a sand tract parallel to the one shown in Figure 8, 
situated further inland (i.e. closer to the sand-trapping region), this asymmetry 
is seen to have become extremely pronounced: the left prongs of the horseshoes 
are only rudimentary, hardly discernable at all in some cases, while the right prongs 
now completely dominate the pattern. From the ground these elongated dunes 
look like those shown in Figure 6b. 


Prototype Specimen 6: Superimposed ridges 

The obliquely-aligned elongated dunes are not only a result of diffusive sand traffic 
across the main wind direction, but—as Experiment 6 proved—are also a factor 
contributing in turn to cross-wind sand traffic. This means that the obliquely-aligned 
elongated dunes in the fringe area of the coastal dune belt translate themselves 
inland by a self-accelerating process, and eventually outrun the rate of sand supply 
coming cross-wind from the coast. Thus the outlying portion of the coastal dune belt, 
if it consists of obliquely aligned dunes, is likely to sever itself at some stage from 
the main body of coastal sand, and then develop a morphology suiting the new 
conditions. Such a morphology is shown in Figure 10, a vertical aerial photograph of 
a sand tract parallel and situated landward to the areas pictured in Figures 8 and 9 
but not adjacent to them: there is a vegetated area of several kilometers width 
between. There is no doubt that this tract was once continuous with the tracts 
described before. Then it developed its basic morphology—and thereafter was severed 
by its own action or by an accident of nature. 

The repetitive direction of the small ridges that are seen superimposed on the 
elongated dune’s ridges and troughs follows the alignment of the pre-existing elon- 
gated dunes. This fact indicates that sand transport takes place not only in the direction 
of the main wind, but is deflected in the landward direction, obviously because of 
the asymmetry of the pre-existing morphology. At the same time it is seen—from 
the alignment of the individual faces of these same ridges, which is strictly transversal 
to the dominant wind direction—that there is a tendency in the system—starting 
with the smaller features—to revert to a symmetrical pattern since the external 
factors of asymmetry (left-hand side source) have ceased to exist. 

This is in line with the observations of Experiment Series 6. 
~ * In some publications these dunes have been referred to as “Seif’’ dunes, while horseshoe dunes 


have been called ‘‘Barchan’” dunes. The author feels that these terms should remain reserved for 
desert type dunes, and not be used for coastal type dunes, in order to avoid confusion. 
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THE HELP OF DYNAMICAL MODELS FOR PURPOSES OF GEOMORPHOLOGICAL 
ANEMOMETRY AND PALAEO-ANEMOMETRY 


One of the purposes of geomorphological studies of sand tracts is to obtain a scientific 
tool for air-photographical evaluation of wind regimes of inaccessible desert areas, 
or, from analysis of fossil landscapes, for the evaluation of wind regimes of past 
geological epochs. Rosenan (1943) and Striem (1954), who have attempted to correlate 
present-day wind regime and dune alignment for sand areas near the ones described 
here, have not reached entirely identical conclusions. The wind regime of past periods 
was discussed by Picard (1943) on the basis of the presence of a huge cemented sand 
(“Kurkar’”’) ridge presumed by him to be a fossil “Seif” dune. His deductions are, 
unfortunately, subject to re-examination because of a misquotation on which he 
based his reasoning*. These attempts in the past so far demonstrate that the problem 
has some urgency, but to be solved, more factors must be taken into consideration 
than could be done by those authors at that time. It is impossible to correlate merely 
main-wind direction and dune alignment. The inter-relation between these and all the 
additional factors (sand supply, sand kind, pre-existing morphology, and traps, 
wind steadiness, etc.) must be studied first. Such a study is possible only by exhaustive 
experimentation, which in turn is possible only by means of dynamical models. 
The selected experiments described in this paper have revealed some of the significance 
of such factors as sand source or sand trap location in relation to wind direction, and 
of pre-existing morphology. But certainly there are other factors operative, not 
included in the present study (e.g. sand composition, treated by the author (1951). 
This state of affairs makes geomorphological anemometry and palaeo-anemometry 
more complicated. But when finally developed, it will render also information 
regarding the environment, in addition to wind-data. 
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Figure 1 
View of experimental flume 


Figure 2 
Miniature landscape photographed through the 
flume window during swift current flow 
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Figure 3 Figure 4 
Miniature landscape photographed through the Miniature landscape photographed through the 
flume window during medium current flow flume window during slow current flow 
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Figure 6a 


Centre part of horseshoe 
dune 


Figure 6b 


Elongated part of 
dune 
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Figure 9 
Vertical aerial photograph of dune landscape adjacent to the one shown in Figure 8 
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Figure 10 
Vertical aerial photograph of dune landscape further inland in the south of Israel 


TENTATIVE ENERGY AND WATER BALANCES FOR THE DEAD SEA 


J. NEUMANN 


Meteorological Service, Israel Ministry of Transport and Communications* 


ABSTRACT 


The average annual rates of evaporation from the northern and southern basins of the 
Dead Sea are calculated from the energy balance of the sea. The resulting evaporation 
rates are as follows: Northern basin—150 cm yr~! (nominally, 147 cm yr~), southern 
basin — 180 cm yr~!. For the Dead Sea as a whole, the rate works out as 155 cm yr}, 

Next, from a combination of the energy balance and water balance approaches, average 
monthly evaporation rates are derived. These figures and figures for monthly lake level 
change lead to estimates of monthly inflow of water into the lake. Finally, the figures for 
monthly inflow are partitioned into component terms in a tentative manner. One of the 
conclusions of the analysis is that the underground inflow of water to the Dead Sea is 
at most of the order of 107 m3 yr“!. 


1. INTRODUCTION 


1.1 Method of calculation 


The annual rate of evaporation from the Dead Sea is calculated in this study by the 
energy balance method. No details of the method will be given here. A review and 
critical remarks on applicational aspects of the method are found in a paper by 
Anderson (1952). A brief review and discussion are given in a paper by the present 
writer (1953a) in connection with an investigation of the energy balances of and 
evaporation from Lake Hula and Lake Tiberias. By considering the energy balance 
of the Dead Sea, the present study completes the investigation for the lakes of the 
Jordan Rift. 

The fundamental idea of the method is that in the process of evaporation, heat 
(latent heat of vaporisation and in the case of a saline lake such as the Dead Sea, 
“differential heat of solution”) is being removed from the water body. By virtue of 
the principle of conservation of energy, it is possible to express this heat energy by 
means of the other terms constituting the energy or heat balance of the water body. 
The terms of the heat balance may be divided conveniently into two groups: radia- 
tion group and convective-conductive group. Evaporation is by far the most impor- 
tant term of the second group. 

The magnitudes of various radiation terms forming the radiation balance of the 
Dead Sea were computed in almost exactly the same manner as followed in the 


* Present address: Department of Meteorology, The Hebrew University of Jerusalem. 


Received October 18, 1955 
137 Bull. Res. Counc. of Israel, Vol. '7G, 1958 


138 J. NEUMANN Bull. Res. Counc. of Israel} 


study on Lake Tiberias and Lake Hula (Neumann 1953a). The only point o 
difference, a minor one, is in the method of correcting for the effects of clouds on 
atmospheric radiation. In this paper, the corrections follow Bolz’s (1949) results. 

All calculations were carried out separately for the northern and southern basins) 
of the Dead Sea (see Figure 1), using the limnological and meteorological data) 


listed in Table I. The meteorological data are those of stations on the coast or near} 
the coast of the lake. Monthly values of the radiation terms are presented in Table II. } 


1.2 Remarks on data 
1.21 Meteorological data. On inspecting Table I, a potential shortcoming becomes | 
apparent immediately: The data for the two basins are not strictly simultaneous and } 
cover but five years each. However, our views of the situation improve if we note |} 
the restricted variability of averages, and particularly those of yearly averages, for ! 
the Dead Sea region. In support of the foregoing statement, we shall compare some 

of the data in Table I with data available through the work of D. Ashbel of the |f 


TABLE I 
Meteorological and limnological data for the Dead Sea 


Tw = lake surface temperature Ta = air temperature 
@a = water vapour pressure in air P = atmospheric pressure 
B =specific gravity of surface waters C =cloud amount 


€w = vapour pressure at lake surface (corrected for the effect of salinity) 
(a) Northern basin 


Tw** Hips Cwt t** €a* P* G~ es 


Month a "C mb mb mb tenths g*cm73 
January 17.0 15.1 16.7 10.9 1066 5.0 
February 16.5 1 16.2 10.9 1065 4.6 
March 18.5 17.2 18.4 11.5 1063 4.4 
April 2255 21.3 23.5 13.3 1061 3.0 
May 26.0 26.8 29.1 16.8 1057 2.8 
June 29.0 29.5 34.5 18.5 1054 0.6 
July 31.9 Bice 38.7 20.9 1051 0.7 1.17 
August 31.5 31.2 39.9 DN2 1052. 0.4 
September 30.5 29.6 37.5 20.6 1055 0.9 
October 28.0 26.9 32.4 18.0 1060 2.9 
November 24.5 22.4 26.5 16.0 1062 4.3 
December 21.0 skys 21.6 12.6 1065 5.0 
re ce a a a RE ee, ee NY PE ee 
Annual average 24.7 23.6 26.8 15.9 1059 2.9 


* Data for Kallia, 1942/46 
** Interpolation between water surface temperatures for the southern basin of the Dead Sea 
and Lake Tiberias (Ashbel 1945, 1950) 
*** Estimate 
**** The annual average of ew is not the average of monthly values. It is the vapour pressure, | 


corrected for the effect of salinity, corresponding to the annual average of water surface 
temperature 
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(B) Southern basin 


Month cone ds Cyt *** @a* p** Cc* (oh AR 
C C mb mb mb tenths g*cm-3 
January 17.5 15.9 15.2 10.5 1066 3.6 1.188 
February 18.7 16.8 16.8 9.6 1065 38) 1.185 
March 215) 19.4 20.1 10.5 1063 2.4 1.183 
April 24.8 24.0 24.4 11.3 1061 1.6 1.182 
May Zit 29.0 28.4 14.1 1057 2.1 1.182 
June 30.2 31.6 33.0 15.7 1054 0.3 1.183 
July 32.0 33.6 36.8 17.3 1051 0.2 1.184 
August 322 33.8 37.1 18.0 1052 0.2 1.185 
September 30.3 31.8 32.4 19.3 1055 0.5 1.187 
October 27.0 28.0 Qed 16.7 1060 1.3 1.188 
November 23.0 23.0 20.9 14.4 1062 Dp 1.189 
December 18.5 18.1 16.3 11.6 1065 Su 1.19 
Annual average 25.28 25.4 24.6 14.1 1059 1.8 1.185 


* Data for Dead Sea Works (Factory), South, 1943/47 
** Data for Kallia, 1942/46 
*** Data due to Ashbel (1950 and files of data) 
**** See note No. 4 for the northern basin 


Department of Meteorology and Climatology, The Hebrew University of Jerusalem, 
covering a longer period. Temperature data for Dead Sea Works, North, are published 
on p. 13 of the supplement Meteorological Normals of Israel to Ashbel’s (1951) 
Regional Climatology of Israel. Data of cloud amount for the latter station as 
well as the temperature data for Dead Sea Works, South, (““Camp”’), are taken from 
files of Ashbel’s data in the archives of the Meteorological Service. 


a) Air temperatures — yearly averages 


Station Period Air Supervising Remarks 
temperature (C) authority 

Kallia 1942/46 23.6 P. M.S. See Table IA 
Dead Sea Works, D.M.C.H.U.; Station 3.5 km 

North 1928/47 ae from 1942:P.M.S. from Kallia 
Dead Sea Works, 

South 1943/47 25.4 P. M.S. See Table IB 
Dead Sea Works, 

South (Camp) 1936/47 25.0 D. M. C. H. U. 


b) Cloud amount — yearly averages 


Kallia 1942/46 2.9 tenths P. M.S. See Table IA 


Dead Sea Works, North 1936/47 3.0 tenths D. M. C. H. U.; 
from 1942: P. M. S. 


P. M.S. = Palestine Meteorological Service ane 
D. M.C.H.U. = Department of Meteorology and Climatology, The Hebrew University of Jerusalem 
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It is seen from the above comparison that the differences between the averages of | 
the data in Table I and those for longer periods are small. For that reason, we may | 


consider the data in Table I as representative for a period of some length. 
A further point illustrating the relatively small variability of meteorological 
elements and their effects on evaporation is offered by the recent work of C. E. Jacob 


(see Section 5 below) who computed the evaporation rate from the Dead Sea by I! 
studying fluctuations of meteorological and hydrometric data, for the years 1930/37 | 


and 1938/44, respectively. The resulting rates differ by Jess than 1%. 


1.22. Temperature of surface waters. Systematic measurements are available only ] 


for the southern basin, initiated and supervised by Ashbel (1950, p. 76, Figure 136). 
Data for the northern basin, which is the major basin (area ca. 760 km? vs. an area 
of ca. 260 km2 for the southern basin, at a lake level of about —392 m M. S. L; 


see de Leeuw, 1946), are deficient. The figures adopted (Table IA) are the outcome | 


of interpolation between surface water temperature data for the southern basin, on 


the one hand, and for Lake Tiberias (Ashbel 1945) on the other hand. In making | 


the interpolation, it was appreciated that deep water bodies (the northern basin has 
an average depth of 185 m, the southern basin 6.4 m; see de Leeuw 1946) tend to 
have somewhat high temperatures in autumn and winter, and somewhat low tempe- 
ratures in spring and summer. Also, account was taken of such scattered data as 
are on hand for the northern basin. 


The annual averages work out as follows: 


Lake Tiberias 23.4°C 
Dead Sea, northern basin 24.7°C 
Dead Sea, southern basin 20.2 C 


1.23 Specific gravity of surface waters. Systematic measurements are available for 
the southern basin (Ashbel 1951, p. 168 and p. 181; 1950, p. 78). Again, for the 
more important northern basin, the situation is far from satisfactory. Some values, 
especially those obtained in the vicinity of the Dead Sea Works, North (see Figure 1), 
and in the western half of the basin, are high (verbal communication by F. Yaron), 
such as 1.18 and even 1.19 g*cm-3, and several of the occasional observations taken 
by Ashbel (1950, p. 77) are similarly high. However, some sporadic measurements 
by Ashbel and by others (Ashbel 1951, p. 180; 1950, p. 11) are low: 1.14 g*cm-3 
and even lower. 


The view has been expressed in discussions that as a consequence of inflow of | 
large volumes of sweet water, the average specific gravity of the surface waters in the | 
northern basin is likely to be appreciably less than in the southern basin. Paucity 
of deta makes it hard to decide the issue, but the following considerations tend to 


disfavour the view just cited. 
The rivers and streams draining to the Dead Sea may be regarded to be like narrow | 


jets of fluid entering a fluid of high density at rest. Now, from an examination of the | 
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Figure 1 
The lower Jordan and the Dead Sea 
D.S. W. = Dead Sea Works P = Police Station (near Kallia) 
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theoretical literature on jets (e.g. Goldstein 1938, Vol. 2, Art. 255), the writer is 
inclined to believe that the area of waters of low specific gravity near the affluence 
of rivers and streams is small. The literature on jets of a temperature and 
velocity different from those of the resting fluid would indicate that the region of 
disturbance in the resting fluid is limited in extent as a result of the mixing deve- 
loping about the edges of jets. 

We have assumed that a representative figure for the average annual specific 
gravity of the surface waters in the northern basin is 1.17 g*cm-3. For the southern 
basin, Ashbel’s figures (see Table I) give an annual average of 1.185 g*cm-3 which 
we adopt. 

As the above figures, particularly the one for the northern waters, may be 
slightly in error, it is worth pointing out that a+ 1% change in the value of specific 
gravity entails a change of + 1%, approximately, in the evaporation rate. Raoult’s law 


Niele Bric es, 2 (1) 
€wo 
where €wo = vapour pressure of pure solvent, e~ = vapour pressure of solution at 
the same temperature and & = mole fraction of solute, does, of course, not hold 
for brines of high concentration, such as the waters of the Dead Sea. However, 
inspection of a recent synopsis by Harbeck (1955, Figure 1) of data on vapour-—pres- 
sure lowering by salts, taken mainly from the “International Critical Tables,’ shows 
that one may write for a not too narrow range of mole fractions, say &’ << &”, 


VA 
ew 


ates 
Sms a SE), (2) 


Cw 


€w= being the vapour pressure of the solution when the mole fraction of the solute 
is &’ and c is a constant of proportionality. From (2), if 4 denotes a finite, but not 
too great change, 


— Ae oc Aé, (3) 


i.e. the vapour pressure of the solution changes /inearly with the concentration é of the 
solute in the range &’ < § < &”. On the other hand, the evaporation rate E can 
be represented by the (empirical) equation 


E= Kaléw — ea), (4) 


Ka being the “coefficient of vapour transfer’’ for the height at which eg, the water- 
vapour pressure in air, is measured. It is clear from equations (3) and (4) that, for 
fixed values of Kg and ea, the variation of E is linear, to a close approximation, with 
ae fraction of the solute or, what comes to the same, with specific gravity of the 
solution. 
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We infer from the above discussion that the evaporation rates obtained below are 
not likely to be in significant error on account of a minor uncertainty in the value 
assigned to the average specific gravity. 


1.24 Vapour pressure of surface waters of the Dead Sea. For an average annual 
cycle, changes in storage of heat in the lake should vanish. Moreover, for an extensive 
water surface, horizontal advection of heat may be considered negligible. Under 
such conditions, the convective-conductive transfer of heat to the atmosphere, 
say Q,, will equal the net radiation surplus Q, at the lake surface, or 


QO, = Q,. (5) 
But 


QO, = 0e-+ Qe, (6) 


Q- standing for the heat utilised in evaporation, which in the case of brines includes, 
as pointed out earlier, in addition to the latent heat of vaporisation and a small 
term representing “energy advected by the evaporated water’’ (for the latter, see 
Anderson 1952), the differential heat of solution (e.g. Glasstone 1951, p. 241). In 
(6), Qc is the symbol for the “‘sensible”’ heat transferred from the lake surface to the 
air. Further, Bowen’s (Bowen 1926, Anderson 1952) ratio is 


Oe _R (7) 
Oc 
LT = Te 
RSeae (8) 
Cw — Ca 


where P = atmospheric pressure in mb, TJ, and T, = temperature of water surface 
and air, respectively, in °C, ew = vapour pressure of the water surface and eg = water 
vapour pressure in air at the height at which Tz is measured; ew» and ég are in mb. 

It is seen from equations (5) to (8) that for partitioning Q, = Q, into its constituents 
QO. and Q¢, we have to know ew, the vapour pressure at the water surface. In the 
case of brines, ey» is, of course, lower than the vapour pressure at the surface of 
distilled water of the same temperature. From the point of view of turbulence theory, 
it is this lowered vapour pressure which accounts for the reduced evaporation rates 
from brines.* 

A group of investigators (M. R. Bloch, Y. Schnerb and F. Yaron), representing the 
Research Council of Israel and the Dead Sea Works Research Laboratory, measured 
the vapour pressure of Dead Sea brine for a range of specific gravities and tempera- 
tures. The results of this work were kindly made available for the present study. 


* On heat-balance theory, the reduced evaporation rate from brines is accounted for by the need 
of supplying the evaporating molecules with the differential heat of solution in addition to the latent 


heat of vaporisation. 
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Most of these vapour-pressure determinations were made for specific gravities and 
temperatures found in the “solar evaporation pans’’ or “salt pans’’ (where the brine 
is evaporated for the purpose of production of salts), which are much higher than those 
occurring in the surface waters of the Dead Sea. In view of this situation, we have 
decided to check how well the experimental data which correspond to natural 
conditions at the surface compare with data for (concentrated) sea water at similar 
chlorinities and temperatures. 


Data of vapour pressure of sea water, for temperatures upward of 25°C and for 
concentrations up to saturation, are available through the work of the Japanese 
chemists Higashi, Nekamura and Hara (see Arons and Kientzler 1954). Arons and 
Kientzler* (1954) have lately completed and perfected the work of Higashi et al. 
by measuring the vapour pressure of sea water from 30°C downward to freezing 
point. In these experiments, chlorinities of sea water varied from 1% to 16%. 
Inter alia, the two authors present a diagram (their Figure 5) which combines, after 
some smoothing, their own results with those of the above-named Japanese chemists. 


On plotting the data obtained for Dead Sea brine by the Research Council-Dead 
Sea Works Research Laboratory group in the diagram representing the smoothed 
results of Higashi et al. and those of Arons and Kientzler for sea water, no systematic 
separation of the data was apparent. It was for this reason that we have decided 
to proceed on the assumption that the lowering of vapour pressure of Dead Sea 
surface waters is the same as for sea water of similar chlorinity and temperature. 
For the northern basin of the Dead Sea, we have taken a chlorinity of 13.5°% (~spe- 


cific gravity 1.17 g*cm-3) and for the southern basin, a chlorinity of 15% (~specific 
gravity 1.185 g*cm-3). 


2. AVERAGE ANNUAL RATES OF EVAPORATION 


According to the results listed in Table II, the average annual radiation balances of 
the two basins are as follows: 


Radiation terms Northern basin Southern basin 


(gcal cm~2 yr—4) (gcal cm? yr—) 

Insolation + 186,000 + 202,000 
Insolation reflected from lake surface — 11,000 — 12,000 
Atmospheric radiation + 276,000 + 281,000 
Atmospheric radiation reflected from 

lake surface — 8,000 — 9,000 
Thermal radiation from lake surface — 333,000 — 335,000 
Net radiation surplus at lake surface + 110,000 + 127,000 


_ ye 


* The writer is indebted to Dr. A. H. Woodcock of Woods Hole Oceanographic Institution, 
Woods Hole, Mass., for making a copy of a report by Arons and Kientzler available in advance of 
publication (Wood Hole Oceanographic Institution, Reference No. 52-94, dated October 1952). 


pe 
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It is seen from the above that the amount of insolation reaching the northern 
basin is estimated at 186,000 gcal cm-2 yr-!. For the southern basin, the corresponding 
figure comes to 202,000 gcal cm-2 yr-!. The reason for this 8% difference resides 
first and foremostin the difference between cloud coverage for the two basins. Table I 
indicates that the average annual cloud amount at Kallia, in the north (see Figure 1), 
is 2.9 tenths, whereas at Dead Sea Works, South (Factory), in the south (see Figure 1), 
the parallel figure is 1.8. A recent analysis (Neumann 1954a) of the relationship 
between insolation and cloud amount shows that, under average conditions, each 
tenth of the sky obscured by clouds is responsible for a reduction by 5 to 6% of 
the insolation which would reach the surface from a cloudless sky.* Estimating the 
cloud coverage over the southern basin as just over one-tenth less compared with 
the northern basin, the calculations thus lead to an insolation rate ca. 6% higher 
than for the north. The remaining 2% of excess in insolation is explained as the 
result of two factors. The first is the difference in solar radiation reaching the top 
of the atmosphere over the two basins. It may be shown, by readily available equations 
or by reference to tables, that the outer atmosphere over the southern basin receives 
some 1,500 gcal cm-? yr-! more solar radiation than over the northern basin. The 
second factor is the lower water-vapour content of the atmosphere in the southern 
region of the lake. A dry atmosphere depletes by absorption less of the solar radiation 
than a humid atmosphere. Table I shows that at the southern station, the absolute 
humidity is 14 mb, while at the northern station it amounts to 16 mb. These are 
surface data, but under average conditions, there is a close connection between the 
water-vapour content of the air near the surface and the water-vapour content of 
the overlying atmospheric column. 


The long-wave radiation balances (= atmospheric radiation Jess atmospheric 
radiation reflected from the lake surface Jess “thermal” radiation emitted by the lake 
surface) of the two basins are almost equal. A consequence of this is that the 
difference in net radiation surplus for the two basins practically equals the difference 
in insolation. 


The net radiation surplus Q,; works out as 110,000 gcal cm-2 yr-! for the northern 
basin and 127,000 gcal cm-2 yr-! for the southern basin. Using equation (8) and 
the data of annual averages in Table I, the Bowen’s ratios R = Q, / Q- [equation (7)], 


and the terms Q- and Qe, Qc + Qe= Q, [equations (5) and (6)], turn out as follows 
(rounded off): 


* By far the major part of the reduction is caused by the reflection of solar radiation from cloud 
tops back to space. Currently, the best estimate (of the world-wide average) of cloud albedo is 0.55 
(private communication from Dr. S. Fritz, U.S. Weather Bureau, Washington, D.C.; see also Fritz 
(1949), Houghton (1954) and Neumann (1954b)). That is to say, 55% of the insolation incident on 
cloud tops is reflected from there. We may add here that the fictive quantity ‘insolation from a 


clear sky” which enters the radiation calculations is ca. 225,000 gcal cm72 yr! for the Dead 
Sea area. 
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Northern basin Southern basin 
R + 0.07 — 0.01 
Qc, gcal cm~2 yr! + 7,000 — 1,000 (10) 
Qe gcal cm~2 yr71 + 103,000 + 128,000 


A negative Bowen’s ratio means that the net transfer of sensible heat is from the 
warmer air to the lake. This heat adds to the net radiation surplus and increases 
the rate of evaporation. 

We have in (10) estimates of the quantity of heat expended by each of the two 
basins of the Dead Sea on evaporation. To express the rate of evaporation in terms 
of depth of water layer we have to remember that in getting 1 g of water vapour by 
evaporation from a brine, we have to supply the brine, as was pointed out before, 
with the so-called ‘differential heat of solution,” in addition to the latent heat of 
vaporisation. Also, there is a minor correction to be made for the heat advected out 
by the evaporated water (for particulars, see Anderson 1952). The latter term we 
estimate at ca. 10 gcal per g of evaporated water (see Neumann 1953a). In brief, the 
evaporation rate E, in cm, is (Anderson 1952, Neumann 1953a) 


Oc 


niet oy 


where 0 is the density of brine in g cm~3 and L is the latent heat of vaporisation, 
in geal g-1, for distilled water at the same temperature. For specific gravities of 
1.17 and 1.185 g*cm~-3, the densities respectively are 1.195 and 1.21 g cm-3. For 
a water temperature of 25°C, such as is the case for the Dead Sea, L = 583 gcal g-!. 
Inserting the appropriate data into equation (11), we arrive at the following 
evaporation rates (N = northern basin, S = southern basin): 


Ey = 147 cmyr“!, E, = 180 cm yr-!. (12) 


To obtain the evaporation rate from the lake as a whole, we have to give a 
weight of 3 to the figure for the northern basin, since its surface area is three 
times greater than the surface area of the southern basin, or 


Exead Sea = 155 em yr}. (13) 


This figure will be used throughout the subsequent sections. 


3. AVERAGE MONTHLY RATES OF EVAPORATION 


The heat of the net radiation surplus at a water surface is disposed in two principal 
ways: part of it is conducted to the atmosphere and the remaining part is conducted 
into the subsurface layers of water where it enters into transitory storage. Unless 
depth-temperature data are available for the water body, itis not generally possible 


148 J. NEUMANN Bull. Res. Counc. of Israel 


to say how much of the radiation surplus was diverted to evaporation and sensible 
heat transfer to the atmosphere. The heat entering into storage, say Q,, may, be 
represented by the equation 

eT w 

Ot 


0, = Co (14) 
where c is specific heat, 9 = density, T» = average temperature of a water column 
and ¢t = time. Clearly, when Ty is a minimum or a maximum, OQ, = 0. This occurs 
twice in the course of a smoothed annual cycle. Once, when the lake reaches its 
lowest temperature about the end of winter, and again, when the lake attains its 
highest ternperature about the end of summer. 

We assume for the Dead Sea that the two months when QO, = 0, are February 
and August, respectively. In adopting these months we are guided by Ashbel’s 
(1945, see also Neumann 1953a) data for Lake Tiberias, which indicate that during 
February and August QO, is practically zero, and by the view that lakes situated in 
one and the same type of climate (and this is approximately true for Lake Tiberias 
and the Dead Sea), tend to have the same thermal behaviour. 

By computing the evaporation rates for February and August, we shall be able 
to compute the rates for the remaining months. In this, we shall utilise the figure for 
the annual rate and some recent results concerning heat-balance theory of water 
bodies. The computation of monthly rates will be completed in five steps. 


3.1 Step No. 1: February and August 


Assuming Q, = 0, then taking values of Q; from Table II and computing the Bowen 
ratios (equation (8)) for February and August, we arrive at the following figures: 


Northern basin Southern basin Weighted average 
(cm) (cm) (cm) 
(rounded off) 
sul sig aie (15) 
February 5.4 6.6 6.0 
August 20.0 22.0 20.0 


3.2 Step No. 2: June, July, September and October 


The Dead Sea and its catchment are situated in a region where no noteworthy wide- 
spread precipitation falls between May and November. A few weeks after the last 
early spring storms, there is no storm-water inflow into the sea. The streams and 
springs settle down to a flow rate which remains practically constant throughout 
summer and early autumn. In hydrology, this flow rate is sometimes designated 
“dry-weather flow” or “base flow”. As examples, we cite figures for the Yarmuk 
stream near its confluence with the Jordan and for the Jordan at Allenby Bridge 
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(see Figure 1). The data for the Yarmuk are from Ionides (1939, p. 139)*, those for 
the Jordan were obtained through the courtesy of the Hydrological Service, Israel 
Ministry of Agriculture. 


TAELE III 


Monthly discharge rates, in m3, of the Yarmuk near its confluence with the Jordan, and the Jordan at 
Allenby Bridge 


2 een ap a ae ae NC NS Seen a 


Jordan 
Yarmuk at Allenby Bridge 
May 1926 — Sep. 1938 Oct. 1935 — Sep. 1947 
(m3 x 10°) (m3 x 109) 
January 32 131 
February 151 143 
March 48 116 
April 34 83 
May 29 84 
June 22 77 
July 21 719 
August 20 83 
September 20 80 
October 21 87 
November 23 87 
December 29 80 
Annual 470 1,135 


It may be objected that the constancy of flow in the Jordan in summer, as indicated 
by the data in Table III, was a result of regulation of this river at Degania by the 
Palestine Electric Corporation. This objection is not quite valid because the main 
season of flow control was in winter (Jonides 1939, p. 137). 

We shall assume that from June to October inclusive the total inflow of water, 
surface and underground, into the Dead Sea takes place at a constant monthly rate J. 
But, for a lake without an outlet 


I=E+AL, (16) 
E being the evaporation rate and AL the amount of change in lake level, both for 


the same monthly period to which J refers. Writing a subscript to indicate the serial 
of the month, our assumption is that Ig = J, = Ig = Ig = Jo, and from (16): 


Eg + ALg = Ey + AL7 = Eg + ALg = Ey + ALy = Ej0 + ALt0 (17) 


Given data for lake level change and the evaporation rate for one of the five months 
- in question, we can at once obtain the evaporation rate for each of the remaining 


four months. 


* In his report on the water resources of Transjordan, Ionides (1939, p. 189) uses the term “‘ground- 
water flow” for what we call “dry-weather flow”, or “base flow”. 
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An estimate of the evaporation rate in August was deduced from heat balance 


considerations and is stated under (15). Data of lake level are also available. Ashbel || 


(1951) published monthly data for the years 1929 to 1946 (p. 181 in his book) and 
a diagram (p. 168) showing the average annual variation. The relevant observations 
were taken near the ‘““Camp”’ of Dead Sea Works, South (Factory). 

For the northern basin, the Hydrological Service (1950, p. 230) published monthly 
data for the years 1939 to 1947. The relevant observations were conducted by the 


Palestine Potash Ltd. in the neighbourhood of their former works in the north (Pales- | 


tine Grid E 199.00, N 131.30). At this location, the gauge was protected from waves. 

There are also two series of unpublished data, copies of records of which are 
deposited in the archives of the Meteorological Service. The first series is by the 
Palestine Potash Ltd., taken near their works in the north. Records of the observations 
go back to 1930. The second series is due to the Public Works Department of the 


former Palestine Government. Records of the data cover the years 1935 to 1947. | 


The place of observations was the former Palestine Police Station situated near the 
Dead Sea Works, North. 

The four sets of data are in qualitative agreement. All four give the highest lake 
level as occurring about the end of April or beginning of May, and the lowest level 
about the middle of November or beginning of December. From the four sets we 
*‘synthesised’’ an annual curve of lake level (Figure 2), according greater weight to 
the data for the north. Special attention was paid to data for the years 1937 to 1944, 
during which the level of the Dead Sea remained constant (see, e.g. Ashbel 1950, 


p. 78, Figure 141), apart from seasonal oscillations. The resulting monthly changes. |] 


in lake level are listed in Table IV. 


TABLE IV 
Monthly changes of lake level, Dead Sea 
(Based on data by Ashbel, the Palestine Potash Ltd. and the Public Works Department of the 
former Palestine Government) 


Month Level change 
(cm) 
January +12.5 
February +16 
March +11 
April + 3.5 
May — 6.5 
June — 5 
July — 9 
August —10.5 
September — 9.5 
October — 7.5 
November — 3.5 
December 4+ 2.5 
Sa 
Annual 0.0 
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Figure 2 
Annual curve of lake level of Dead Sea 


By utilising the figure for the evaporation rate in August given under (15) (20 cm, 
Eg = 20) and the data in Table IV, equation (16) results in the following evaporation 
rates for the months June, July, September and October: 


E«6 = 14.5 cm 

E; = 18.5 

Ey = 19 th) 
E\o = 17 


3.3 Step No. 3: November 


Consider the idealised case where the annual variations of all three terms of the 
water-balance equation (16) are perfect sine curves. In this case, in the month when 
AL =0, I= E, with 7 and E being respectively the average monthly inflow and 
evaporation, i.e. one-twelfth part of the average annual figures. 

It may be inferred from an energy-balance investigation (Neumann 1954a) of the 
problem of annual variation of evaporation from lakes that this variation comes 
nearest to a sine curve in the autumn and winter seasons. Now, during the autumn- 
winter semester of the year, for the Dead Sea, AL = 0 about the end of November. 
Hence, following what was said in the foregoing paragraph, a good estimate of the 
evaporation rate for the monthly period centred about the end of November would 
_ be one-twelfth part of the estimated annual evaporation rate, i.e. (155/12) = ~ 13 cm. 
But since we are interested in the evaporation rate for November proper, we have 
to increase slightly the above figure. By considering the rate deduced for October, 
we suggest to adopt for November the figure of 13.5 cm, or 


Ey = 13.5 cm. (19) 


3.4 Step No. 4: December and January 


Let us here introduce the term “annual heat budget” Q,,ofa water bodysuchasa lake. 
The annual heat budget is the amount of heat required to raise the temperature 
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of a column of unit cross-section of the lake from its minimum temperature (late) |} 
in winter to its maximum temperature (late) in summer. This heat amount, which |} 
enters into storage in the lake during the spring and summer Seasons, leaves the lake | 
during the autumn-winter semester of an annual cycle. It will be seen from Figure 4 | 
that deep lakes possess an annual heat budget of between 20,000 and 40,000 gcal | 
cm-2, Thus the heat released from storage in lakes of some depth during the | 
cool season is an important source of heat for the lake surface to transfer to the ] 
atmosphere. 

The writer has shown elsewhere (1953b) that Q,, may be computed from surface | 
data alone, without soundings of depth temperatures, provided that an assumption is_ | 
made in regard to the dates when the sea or lake reaches its minimum temperature | | 
and maximum temperature, respectively. Earlier, in “Step No. 1”, we have assumed 
that in the case of the Dcad Sea, these dates occur in February and August. We 
shall now make the more exacting assumption that the dates occur at mid-February | 
and mid-August, respectively. 

As before, we shall denote by Q; the net radiation surplus at the lake surface, 
R = Bowen’s ratio (R = Q-/Q-), L = latent heat of vaporisation, 9 = density of 
surface waters, a = a correction term accounting for the heat advected out by the 
evaporated water (earlier we have assumed for the Dead Sea a = 10 gcal g-!) and 
E = evaporation rate. 

The above symbols marked by a prime will represent the values of these quantities 
for the period between mid-February and mid-August; a double prime will refer 
to the period between mid-August and mid- mpage Without primes, the symbols 
relate to a complete annual cycle. 


It is shown in the paper referred to above, that 


Cite ee ey ee 


and 
: QO, — Oy _ Of’ + On 
UD a a Bt ee a Eee 
elt Rte] vm atR) tay ev 
Naturally, 
E' = Wo a E. (22) 


We shall now apply these equations to the problem of the Dead Sea. Table V 
lists for the northern basin, on the basis of Table I, half-yearly averages of the meteo- 
rological and limnological data which are required for obtaining the half-yearly 
Bowen’s ratios R’ and R”. Also added are the half-yearly net radiation balances Or’ 
and Q,'”. The latter were compiled from Table ITA. 
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TABLE V 
Half-yearly (15 February to 15 August and 16 August to 14 February) averages of meteorological 
and limnological data for the northern basin of the Dead Sea 
i 


= 25.2°C Tw =" 24:2°C 
To = 249°C Tat =" 223°C 
ew = 28.7mb Co A = 1627. limb 
és = 16.2mb Ga’) = gel5:6mb 
hem = 1,057 mb Ret = 1,063 mb 
QQ,’ = 76,000 gcal cm-2 Q,’’ = 34,00 gcal cm—2 
R= .0:02) Re = 0.11 
eo == 119gcem3 Q”’ = 1.19gcm-3 
Eo ==" 583 gcalipy Te = 583 gcalg 

R = 0.07 


ES SE Ses 1s ioe Ere ER ia case (Te See Rea met 
Introducing the data of Table V into equations (20) and (21), it is found that 


QO = 30,000 gcal cm-2 
Ey = 65cm (23) 
B= 82.001, 


N being the symbol for the northern basin. Ey’ + Ey’ = Ey = 147 cm, as ex- 
pected. We note in passing that according to the above results, for the northern 
basin, the evaporation rate in the autumn-winter semester exceeds the evaporation 
rate in the spring-summer semester. This is reasonable in the light of what we know 
about the annual course of evaporation from deep water bodies and in view of 
the great depth of the northern basin (average depth ca. 185m, de Leeuw 1946). We 
shall see further on that the reverse relationship is obtained for the shallow 
southern basin. 

Before being able to utilise the figures under (23), we have to find the corresponding 
figures for the southern basin. It is, however, not possible to follow the procedure 
just described, for in the case of the southern basin, the Bowen’s ratios R’, R’’ and R 
are all small and nearly equal. Their use in conjunction with equations (20) and (21) 
might lead to serious errors. We therefore have to proceed in a different manner. 

Let T be the annual maximum temperature, T, the annual minimum temperature, 
both in °C, o the density in g cm-3 and c the specific heat, in geal g-! deg-!, of the 
water column whose annual heat budget is Q,. Then, we may define an equivalent 
lake-depth D, in cm, such that ' 


D(Tz —Tn) 0°¢ 


on = Re iy tee (24) 


The factor 9 :c may be taken as unity, which is a fair approximation to truth, so 
that, numerically, 


D(Tz — Tn) 
p gerest Maa dla ) 
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In a deep lake, there will be an epilimnion. Most of the heat represented by Oy 
entering into storage in the lake will be confined to the epilimnion. If the lake is 
very deep, as is the case with the northern basin of the Dead Sea, the equivalent 
depth D will be less than lake depth. Nor will it be equal to the depth of the epilimnion. 
In equation (25) the term Q,, may be represented (see Figure 3) as the area of a 
triangle whose base is (Tz — Tn) and its height D. 


AMPLITUDE OF ANNUAL TEMPERATURE WAVE 
SURFACE WATERS (°C) 


- 20000 


EQUIVALENT DEPTH OF BASIN (m9) 
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Figure 3 
Calculation of the annual heat budget of the southern basin from the annual heat budget of the 
northern basin. Triangle PSR = northern basin, trapezoid PP’R’R = southern basin 


Now, if in lakes the temperature range varied with depth linearly, D would be 
equal to epilimnion depth (assuming that heat-content variation is confined to 
the epilimnion). In reality, below an upper homothermal layer, this variation is 
something like an exponential reduction with depth, so that the reduction with 
depth takes place at a rate faster than linear. As a consequence, D should be less, 
probably much less, than epilimnion depth. We have deduced earlier for the northern 
basin of the Dead Sea an annual heat budget Q,, = 30,000 gcal cm-2; the term 
(Tz — Tn) may be obtained, to an approximation, from Table I(A). Taking T; = 
315°C, Ty ='16.5°C; (Tz: — Tn) = ASC. 

Inserting these data* into equation (25), we get D = 4,000 cm or 40 m, whereas 
the depth of epilimnion in the northern basin appears to be ca. 75 m. The latter we 
infer from salinity data due to Blake (1930, p. 7) and information due to Novomeisky 
cited by Braslavsky (1950, p. 15). Salinity data show little variation below about 75 m. 

To estimate Q,, for the southern basin, we return to Figure 3 and equation (25). 
For the southern basin, from Table I(B), Tz = 32.2°C, Tn = 17.5°C, (Tz — oe 


* Strictly speaking 7, = 31.5°C is not the maximum temperature and Tn = 16.5°C is not the 
minimum temperature required. Rather they are average temperatures of the warmest and coldest 
months, respectively. We use them as the proper maximum and minimum temperatures are not 
available. The error incurred is compensated for by the value of D calculated from equation (25). 
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14.7°C. Since this temperature range is nearly the same as that for the northern 
basin, i.e. 15°C, we shall adopt the latter figure. This means that in Figure 3 we can 
use the same “base”’ PR both for the northern and southern basins. But since the 
southern basin is shallow, its average depth being 6.4 m (de Leeuw 1946), the tempera- 
ture range at its bottom will not reduce to zero. Hence, it will be more appropriate 
to represent in Figure 3 the annual heat budget of the southern basin by a trape- 
zoid whose base is PR. 

Coming to the question of height (or, rather, depth) of the trapezoid we note that 
in the case of a shallow lake, allowance must be made for the seasonal variation of 
heat content in the lake bed. Heat enters the bed in spring and summer and, under 
average conditions, the same amount of heat should leave the bed in the remaining 
seasons of a full annual cycle. The alternating part of the lake bed heat content 
may be written 


f5co Te dz, (26) 


where c = specific heat, @ = density, 7, = alternating part of the temperature wave 
at depth z in the lake bed. By following Bowden’s (1948) evaluation of that term 
for the Irish Sea, or by repeating the more detailed considerations by the present 
writer (1953a) concerning Lake Hula, we find that the heat budget of the lake bed 
is equivalent to the heat budget of a layer of water ca. 1 m thick at the bottom of the 
lake. Hence, in order to allow for the effect of lake bed, we should increase the 
figure for lake depth by 1 m. We thus assume, for the purposes of evaluation of Oy 
by means of Figure 3, that D = 740 cm (= 640 + 100). 

It is readily confirmed that for D = 740, the trapezoid PP’R’R in Figure 3 has an 
‘area’, that is to say, an annual heat budget, of ca. 10,500 gcal cm-2, or about 
one-third of the annual heat budget computed for the northern basin. 

A confirmation of the approximate correctness of the above result is offered by 
Figure 4, which plots data of annual heat budgets of over 100 lakes vs. average lake 
depth. Most of the data were obtained through the courtesy of G. E. Hutchinson, 
Osborn Zoological Laboratory, Yale University, New Haven, Conn. The remaining 
data were collected by the writer from limnological literature on hand. The free-hand 
curve drawn in Figure 4 represents a “‘best’’ estimate of annual heat budget from 
data of average depth for lakes. It will be seen that for an average depth of 6.4 m, 
the diagram gives an annual heat budget Q;, = 10,000 gcal cm~2. 

We adopt for the annual heat budget of the southern basin of the Dead Sea the 
round figure of 10,000 gcal cm-2 and use it in computing Ey’, where the subscript 
symbolises the southern basin. By having estimates of E” both for the northern and 
southern basins, we shall be able to compute E” for the lake as a whole and arrive 
at estimates of the evaporation rates in December and senuety ek in spring — see 


Sub-section 3.5). 
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Figure 4 


The annual heat budget of lakes as a function of average depth of the lake. (Most of the data was 

obtained through the courtesy of G.E. Hutchinson, Osborn Zoological Laboratory, Yale University, 

New Haven, Conn.) The names of some of the better-known lakes are stated in the lower right-hand 
corner 


The data needed for computing E,” are listed in Table VI. They were obtained 
from Table I(B). 


With the data on the right hand side of Table VI, equation (21) gives 


E.. = 130in, (27) 
and since for the southern basin, E = 180 cm, 
E, = 107 cm. (28) 


Thus in the southern basin in spring and summer, the evaporation rate considerably 
exceeds the rate for the autumn-winter semester. This is the opposite of what we 


TABLE VI 
Meteorological and limnological data for the southern basin of the Dead Sea for the period from mid- 
August to mid-February 
ce ee ee, AO AE ONS TIES ECE ONE 


Fe” = "BAC r’ = 42,000 geal cm~2 
Ts’ = 26.6°C R”’ = —0.01 

ew’ = 26.3 mb L’ = 585 gcalg7! 
es’ = 14.4mb eo” = 121gem> 

Pp = 1063mb a = 10gcalg™ 


te ee 
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found for the northern basin. Far from being at variance, the two results are in accord 
with what we know about the annual variation of evaporation from shallow and 
deep water bodies, respectively. Once more denoting by subscript N and S the nor- 
thern and southern basins, in that order, we have 


Ey = 147.cm £ == 65'em Eo cnt 29 
E,; = 180cm Be. =O] cm Fao 13m) Ee 

Taking the weighted averages of the E,”, we get for the Dead Sea as whole, 
E” ='80cm. (30) 


So far, we have the evaporation rates for August, September, October, November 
and February. Remembering that E” refers to the period from mid-August to mid- 
February, we shall have to divide the rates for August and February in two parts. 
Of the 20 cm evaporation rate for August (see equation (15)), we tentatively assign 
10.5 cm to the second half of August and of the 6 cm rate for February, we assign 
2.5 cm to the first half of February. Then, in view of equations (18) and (19) we have 


1/2 August 10.5 cm 

- September 19 “cm 
October 17 sem (31) 

November 13.5 cm 

1/2 February 2:5.¢m 


Total 62.5 cm 


For the period as a whole we have, from equation (30), a total evaporation of 80 cm. 
Therefore, the combined amounts of evaporation for December and January should 
come to 17.5 cm. On the basis of graphical interpolation, we suggest to divide that 
quantity as follows: 
ad ' E12 = 10cm 


32 
Ey = 7.5m ed 
3.5 Step No.5: March, April and May 
It follows from (29) that for the Dead Sea as a whole, 
Pve= T5.cit (33) 


Remembering that E’ refers to the period from mid-February to mid-August, 
we have, so far, the following estimates [see (15) and (18)]: 


1/2 February 3.5 cm 
June 14.5 cm 

July 18.5 cm 

naj2 August 9.5 cm 


Total 46.0 cm 
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leaving, in view of (33), a total amount of 29 cm to be divided up among March, | 
April and May. | 

One should not attempt to partition the 29 cm amount by straightforward | 
interpolation. It has been shown elsewhere (Neumann 1954) that evaporation from | 
lakes in spring possesses a peculiar feature. Instead of the evaporation rate increasing, | 
so to speak, monotonically from winter to summer, the initial increase early in spring | 
is followed by a decline in the evaporation rate about the middle of spring or so, to be, | 
in turn, succeeded by a renewed increase. 

The spring minimum of evaporation appears to be a fact for the Dead Sea. All | 
the available level data indicate that the Dead Sea is at its highest level about the 
end of April or early in May (see Ashbel 1951, p. 181 and diagram on p. 168; Hydro- 
logical Service 1950, p. 230; Figure 2 in the present paper). Since the inflow of the 
surface waters reaches its maximum in February (see data for the Jordan in Table III) 
and this maximum is followed by a rapid drop in surface inflow, it is hard to explain — 
the high water in spring by a factor other than a relatively low rate of evaporation 
and, perhaps, by a high rate of underground water inflow. However, the results of 
the present paper do not support conjectures as to the existence of an underground 
inflow in significant amounts (see below), an inference which is in agreement with 
recent work done by C.E. Jacob. We shall return to this point subsequently. 

We shall tentatively partition the evaporation amount of 29 cm as follows: 


March 8.5 cm 
April 9.0 cm (35) 
May 11.5 cm 


but emphasise the uncertainty of these figures. 


3.6 Summary of results 


Collecting results, including the estimates under (35), we arrive at Table VII and 
Figure 5. 
TABLE VII 
Monthly rates of evaporation from the Dead Sea 


Month Cm 
January 7.5 
February 6.0 
March 8.5 
April 9.0 
May 11.5 
June 14.5 
July : 18.5 
August 20.0 
September 19.0 
October 17.0 
November 13.5 
December 10.0 


Annual 155 
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4. TENTATIVE WATER BALANCE FOR THE DEAD SEA 


The results of the foregoing section make it possible to set up a tentative water 
balance for the Dead Sea for an average yearly cycle, assuming no net change in the 
volume of water stored in the lake. The highly tentative character of the attempt 
undertaken in this section is emphasised. 


In Table VIII, column B reproduces figures for monthly changes of lake level 
given in Table IV. It was pointed out earlier that these figures were obtained by 
preparing from all available level data a smoothed, idealised curve. In column A, 
we re-list the monthly rates of evaporation deduced in the present study as summarised 
in Table VII. In view of equation (16), the fundamental water-balance equation 
for a lake without an outlet, the sum of columns A and B represents the total 
inflow. In columns D, E, F and G, we partition column C into the following four 
“components”: discharge of the Jordan at Allenby Bridge (see Figure 1), rainfall 
in the Dead Sea, dry weather flow + underground inflow, and, finally, storm- 
water inflow. 

The figures in column D for the Jordan were prepared on the basis of discharge 
data for this river at Allenby Bridge, kindly made available by the Hydrological 
Service, Israel Ministry of Agriculture. We assume a yearly discharge of 1.25 x 109 m3, 
which for a lake area of 1020 km? (such as is fairly representative for the last 
decades) reduces, in round figures, to 122.5 cm. Since this figure might appear slightly 


TABLE VIII 


Tentative monthly water balances for the Dead Sea 
(all terms in cm) 


A B (S D E F G 
Month Evapora- Change of Total Jordan at Rain- Dry-weather Storm-water 
tion level inflow Allenby fall surface inflow inflow 
Bridge south of 
Allenby Bridge 
and under- 


ground inflow 


January BS) +12.5 20 14.5 2.5 1.25 1.75 
February Gil! +16 22 16 1.5 1.25 3.25 
March 8.5 +11 19.5 14 1 125 3.25 
April 9 + 3.5 12.5 10 0.5 1.25 0.75 
May 11.5 — 0.5 11 9.25 — 1.25 0.5 
June 14.5 — 5 9.5 8.25 _ 1.25 - 
July 18.5 50) 9.5 8.25 ie 1.25 “i 
August 200' —10.5 9.5 8.25 — 1.25 — 
September 19 — 9.5 9.5 8.25 — 1.25 — 
October 17 — 7.5 9.5 8.25 — 1.25 _ 
November 13.5 — 3.5 10 8.25 0.5 1525 _ 
December 10 + 2.5 12.5 9.25 1.5 1.25 0.5 
ul itt gprs Vilmeey ss pap ithe Bea fan hy Eats a 
Annual 155 0 155 122.5 aS 15 10 
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high, it is well to point out that the Jordan’s average yearly discharge in the years 
1937/38 to 1943/44, during which the level of the Dead Sea remained stable, came 
very close to 1.25 x 109 m3. The precipitation estimates in column E are due to 
N. Rosenan. 

The figures in column F are to represent dry-weather surface inflow (or base flow 
at the surface) and underground inflow. We assume that these terms are constant 
throughout the year. It will be remembered that the assumption of constancy of 
these terms enabled us to compute monthly evaporation rates for June, July, Septem- 
ber and October and that in support of the assumption we introduced the data 
in Table III. 

Figures in column G for storm-water inflow are, of course, the balance left over 
from total inflow after subtracting the terms for the Jordan, rainfall and dry weather 
surface flow plus underground inflow. 

Summarising the table, we get, for a lake area of 1020 km2, the following annual 
water balance, in units of 106 m3: 


Gains Losses 


Jordan at Allenby Bridge 1250 Evaporation 1580 
Rainfall 80 
Dry-weather surface inflow south of 

Allenby Bridge and underground 


inflow to Dead Sea 150 
Storm-water inflow 100 
1580 1580 


It is worth adding a remark on the term representing dry-weather surface inflow 
and underground inflow. 

The total surface inflow from the west scarp of the Dead Sea is estimated by 
Blake and Goldschmidt (1947, p. 399; see figures for the Ein Gedi and Ein Fashkha 
groups of springs) at ca. 75 x 106 m3 yr-3. Of this, about 60 x 106 m3 yr-! constitutes 
dry-weather flow, according to an estimate given by M. Goldschmidt at a meeting 
in summer 1954 to which we will refer once more in Section 5. For the eastern coast 
of the lake, a map and a table by Ionides (1939, map preceding p. 147; and his table 
21, p. 190) give a dry-weather flow (“ground water flow”) of nearly 100 x 106 m3 
yr-!. Thus, if we accept the foregoing two figures, we should have a total dry-weather 
surface flow of about 150 x 106 m3 yr-!, leaving very little for an underground inflow. 

Assuming that Ionides’ figures are high, underground inflow should hardly amount 
to more than 50 x 106 m3 yr-1. All these figures are, of course, uncertain to some 


extent, but the conclusion that underground inflow is relatively small seems to be 
inescapable. 
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5. RELATION OF PRESENT RESULTS TO RECENT RESULTS 
BY OTHER INVESTIGATORS 


Ionides (1939, p. 143) states: “It is calculated that the depth evaporated from the 
surface (of the Dead Sea) annually is an amount approaching two metres.” However, 
no details of the calculation are published. A similar evaporation rate is suggested 
by de Leeuw (1946). It is probably correct to say that Ionides, just as de Leeuw, 
was influenced by the assumed existence of a large underground inflow to the Dead 
Sea and if this assumption is relinquished, the estimates both by Ionides and de 
Leeuw would come close to the result reached in the present paper. 


Ashbel (1951, p. 180) estimates that during recent years the total water supply 
to the Dead Sea amounted to between 1.25 and 1.40 x 106 m3 yr-1, or from 125 
to 140 cm yr-!. But since during the very same years, the level of the Dead Sea 
fell, it is clear that the evaporation rate must have exceeded the supply rate. 


The results of the present study bear a close relationship to the results of recent 
work done by C.E. Jacob (Brigham Young University, Provo, Utah) and by M. 
Goldschmidt (Hydrological Service, Israel Ministry of Agriculture). 


The work concerning the annual rates of evaporation from the northern and 
southern basins of the Dead Sea, as summarised in Section 2 above, was completed 
in spring 1954, and reported in June 1954 at a symposium of the Hydrology Section, 

‘Israel Society of Geodesy and Geophysics. About the same time, C.E. Jacob, in 
his capacity as consultant to Water Planning for Israel Ltd., began his work on the 
Dead Sea. Jacob undertook a study of fluctuations in meteorological (air temperature, 
Piche evaporation, rainfall) and hydrometric (Jordan discharge and lake-level data) 
variables. His results, which became available in summer 1954, and which included 
evaporation figures for bi-monthly periods, indicated an evaporation rate of 170 cm 
yr-1, 

Stimulated, in part, by conversations with C.E. Jacob, the present writer extended 
his study to monthly periods. That phase of the study, which is summarised in 
Section 3 above, was completed in August 1954. 


In the course of discussions, Jacob agreed that his figures for the spring months 
should be reduced slightly in order to allow for the phenomenon of “spring minimum 
of evaporation” referred to earlier in Sub-section 3.5. This correction is desirable, 
as the phenomenon of low evaporation rates from lakes in the spring is not reflected 
in air temperature and Piche evaporation data on which his study is based. Jacob's 
corrected result amounts to ca. 160 cm yr-!. 


The result of Jacob’s work as well as the monthly evaporation rates derived by 
the present writer were reported by Jacob at a meeting held in the offices of Water 
Planning for Israel Ltd. at Hakirya in August 1954. At the same meeting, M. Gold- 
schmidt presented his estimates based on hydrometric data and estimates for the 
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ungauged surface inflow. Assuming no underground inflow, he arrives at an 
evaporation rate of 153 cm yr-!. 

Thus the results of Jacob, Goldschmidt and the present writer, developed along 
different lines of approach, are close. It would appear that the average evaporation 
rate is between 150 and 170 cm yr-!. It is quite possible that the rate lies in the narrower 
range of 155 to 160 or 155 to 165 cm yr-1. 
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LETTER TO THE EDITOR 


A note on the revision of the stratigraphy of the Southern Carmel. E. KAsHat, Department of Geology, 
The Hebrew University of Jerusalem. 

Since the publication of our paper ‘“‘On the Lithostratigraphy and tectonics of the Carmel”?, paleon- 
tological investigations have led us to certain stratigraphic revisions discussed briefly in this note and 
illustrated in the attached geological map and section (Figures 1 and 2). 

1. The reef limestone and dolomite beds between the volcanic tuffs appearing on the western 
flank of Wadi Natif east of Umm es Zeinat are covered with chalk referred formerly to the Daliya 
series. The fossils of this chalk determined now, however, indicate an age identical with the Upper 
Cenomanian, so-called Juneidiye, series. 


2. The thin marls around Umm es Zeinat are now considered stratigraphically equal to the thicker 
Turonian Daliya marl. They are no longer referred to the Zikhron Zone. Therefore, the layers under- 
neath the Daliya marls belong to the Muhraqa series considered to be the base of the Turonian. 
The Qumbaze limestone remains the top of the Turonian. Both Turonian horizons, Daliya mar] and 
Muhraqa limestone, are thinning out towards the south so that Daliya marl disappears completely 
towards Ijzim. Muhraqa limestone series is represented from Ijzim to Shune by the reef limestone of 
the lower Shune series. 


3. Therefore the formations underlying the Shune series, as found in most of the Lower and Western 
Carmel up to Wadi Mughrara, are of Cenomanian age and are stratigraphically identical with 
the Main Chalk complex developed principally in the Upper Carmel. Many difficulties arose in the 
comparison of these formations because of the abundant reef facies occurring in the Lower Carmel, 
contrasting with the rather uniform chalky lithology of the Upper Carmel. 


4. These stratigraphic changes naturally lead us also to an age revision of the various volcanic 
intercalations as expressed in the following table. 


2 km north of Bat Shlomo Senonian 
Around Shefeya Shuna series 
Around Ofer Shuna series 
South of Ein Ghazal Shuna series 

1 km south of Ijzim Shuna series 
Around the valley of Maqura Juneidiye series 
The valley east of Ijzim Juneidiye series 
East of Zikhron Ya’akov Juneidiye series 
Umm_et Tos (1 km NE of Fureidis) Juneidiye series 
North of Ein Ghazal Juneidiye series 
Along the road E of Ijzim Juneidiye series 
W and E of Muhraga Khureibe-Juneidiye series 
1 km south of Ein Ghazal Khureibe series 
The valley W of Ijzim Khureibe series 
The valley E of Umm ez Zinat Khureibe series 
Wadi Zarra’a (i km S of Damun) Isfiya series 


eee 
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Figure 1 
Revised stratigraphic map of the southern Carmel. 
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